Changeset 34968 in project
 Timestamp:
 12/29/17 20:58:25 (3 months ago)
 Location:
 release/4/srfi27/trunk
 Files:

 14 added
 3 edited
Legend:
 Unmodified
 Added
 Removed

release/4/srfi27/trunk/srfi27distributions.scm
r34967 r34968 1 1 ;;;; srfi27distributions.scm 2 2 ;;;; Kon Lovett, Dec '17 3 ;;;; Kon Lovett, Jun '174 ;;;; Kon Lovett, May '065 3 6 ; Chicken Generic Arithmetic! (could use fp routines) 7 8 (module srfi27distributions 9 10 (;export 11 makerandomnormals 12 makerandomexponentials 13 makerandomtriangles 14 makerandompoissons 15 makerandombernoullis 16 makerandombinomials 17 makerandomgeometrics 18 makerandomlognormals 19 makerandomcauchys 20 makerandomgammas 21 makerandomerlangs 22 makerandomparetos 23 makerandomlevys 24 makerandomweibulls) 4 (module srfi27distributions () 25 5 26 6 (import scheme chicken) 27 7 8 (reexport 9 srfi27normals 10 srfi27exponentials 11 srfi27triangles 12 srfi27poissons 13 srfi27bernoullis 14 srfi27binomials 15 srfi27geometrics 16 srfi27lognormals 17 srfi27cauchys 18 srfi27gammas 19 srfi27erlangs 20 srfi27paretos 21 srfi27levys 22 srfi27weibulls) 28 23 (use 29 (only typeerrors errorargumenttype) 30 (only typechecks 31 definecheck+errortype 32 checkprocedure 33 checkcardinalinteger 34 checkreal 35 checkopeninterval 36 checkclosedinterval) 37 srfi27 38 srfi27distributionssupport) 39 40 ;;; Normal distribution 41 42 ;; Knuth's "The Art of Computer Programming", Vol. II, 2nd ed., 43 ;; Algorithm P of Section 3.4.1.C. 44 45 (define (*makerandomnormals mu sigma randoms) 46 (let ((next #f)) 47 (lambda () 48 (if next 49 (let ((result next)) 50 (set! next #f) 51 (+ mu (* sigma result))) 52 (let loop () 53 (let* ( 54 (v1 ( (* 2.0 (randoms)) 1.0) ) 55 (v2 ( (* 2.0 (randoms)) 1.0) ) 56 (s (+ (* v1 v1) (* v2 v2)) ) ) 57 ; 58 (if (<= 1.0 s) 59 (loop) 60 (let ((scale (sqrt (/ (* 2.0 (log s)) s)))) 61 (set! next (* scale v2)) 62 (+ mu (* sigma scale v1))))))))) ) 63 64 (define (makerandomnormals #!key (mu 0.0) (sigma 1.0) (randoms (randomreal/current))) 65 (checkreal 'makerandomnormals mu 'mu) 66 (checknonzeroreal 'makerandomnormals sigma 'sigma) 67 (checkprocedure 'makerandomnormals randoms 'randoms) 68 (values 69 (*makerandomnormals mu sigma randoms) 70 (lambda () (values mu sigma randoms))) ) 71 72 ;;; Exponential distribution 73 74 ;; Knuth's "The Art of Computer Programming", Vol. II, 2nd ed., 75 ;; Section 3.4.1.D. 76 77 (define (*makerandomexponentials mu randoms) 78 (if (= 1.0 mu) 79 (lambda () ( (log (randoms)))) 80 (lambda () (* mu ( (log (randoms)))))) ) 81 82 (define (makerandomexponentials #!key (mu 1.0) (randoms (randomreal/current))) 83 (checkrealunit 'makerandomexponentials mu 'mu) 84 (checkprocedure 'makerandomexponentials randoms 'randoms) 85 (values 86 (*makerandomexponentials mu randoms) 87 (lambda () (values mu randoms))) ) 88 89 ;;; Triangle distribution 90 91 ;; s  smallest, m  most probable, l  largest 92 93 (define (*makerandomtriangles s m l randoms) 94 (let ((d1 ( m s)) 95 (d2 ( l s)) 96 (d3 (sqrt ( l m)))) 97 (let ( 98 (q1 (/ d1 d2) ) 99 (p1 (sqrt (* d1 d2)) ) ) 100 ; 101 (lambda () 102 (let ((u (randoms))) 103 (if (<= u q1) 104 (+ s (* p1 (sqrt u))) 105 ( l (* d3 (sqrt ( (* d2 u) d1))))))))) ) 106 107 (define (makerandomtriangles #!key (s 0.0) (m 0.5) (l 1.0) (randoms (randomreal/current))) 108 (checkreal 'makerandomtriangles s 's) 109 (checkreal 'makerandomtriangles m 'm) 110 (checkreal 'makerandomtriangles l 'l) 111 (checkrealopeninterval 'makerandomtriangles l s +inf.0 'l) 112 (checkrealclosedinterval 'makerandomtriangles m s l 'm) 113 (checkprocedure 'makerandomtriangles randoms 'randoms) 114 (values 115 (*makerandomtriangles s m l randoms) 116 (lambda () (values s m l randoms))) ) 117 118 ;;; Poisson distribution 119 120 (define (*makerandompoissons mu randoms) 121 (let ((emu (exp ( mu)))) 122 (lambda () 123 ;FIXME O(mu) but O(log(mu)) desired for >> mu 124 (do ((m 0 (fx+ 1 m)) 125 (prod (randoms) (* prod (randoms)))) 126 ((<= prod emu) m)))) ) 127 128 (define (makerandompoissons #!key (mu 1.0) (randoms (randomreal/current))) 129 (checknonnegativereal 'makerandompoissons mu 'mu) 130 (checkprocedure 'makerandompoissons randoms 'randoms) 131 (values 132 (*makerandompoissons mu randoms) 133 (lambda () (values mu randoms))) ) 134 135 ;;; Bernoulli distribution 136 137 (define (*makerandombernoullis p randoms) 138 (cond 139 ((= 0.0 p) (lambda () #f)) 140 ((= 1.0 p) (lambda () #t)) 141 (else (lambda () (<= (randoms) p)))) ) 142 143 (define (makerandombernoullis #!key (p 0.5) (randoms (randomreal/current))) 144 (checkrealunit 'makerandombernoullis p 'p) 145 (checkprocedure 'makerandombernoullis randoms 'randoms) 146 (values 147 (*makerandombernoullis p randoms) 148 (lambda () (values p randoms))) ) 149 150 ;;; Binomial distribution 151 152 (define (*makerandombinomials t p randoms) 153 (let ((bernoullis (*makerandombernoullis p randoms))) 154 ;FIXME O(t) but O(log(t)) desired for >> t 155 (if (fixnum? t) 156 (lambda () 157 (do ((i 0 (fx+ 1 i)) 158 (n 0 (if (bernoullis) (fx+ 1 n) n))) 159 ((fx<= t i) n))) 160 (lambda () 161 (do ((i 0 (add1 i)) 162 (n 0 (if (bernoullis) (add1 n) n))) 163 ((<= t i) n))))) ) 164 165 (define (makerandombinomials #!key (t 1) (p 0.5) (randoms (randomreal/current))) 166 (checkcardinalinteger 'makerandombinomials t 't) 167 (checkrealunit 'makerandombinomials p 'p) 168 (checkprocedure 'makerandombinomials randoms 'randoms) 169 (values 170 (*makerandombinomials t p randoms) 171 (lambda () (values t p randoms))) ) 172 173 ;;; Geometric distribution 174 175 (define (*makerandomgeometrics p randoms) 176 (let ((logp (log p))) 177 (lambda () 178 (+ 1 (inexact>exact (floor (/ (log ( 1.0 (randoms))) logp)))))) ) 179 180 (define (makerandomgeometrics #!key (p 0.5) (randoms (randomreal/current))) 181 (checkrealunit 'makerandomgeometrics p 'p) 182 (checkprocedure 'makerandomgeometrics randoms 'randoms) 183 (values 184 (*makerandomgeometrics p randoms) 185 (lambda () (values p randoms))) ) 186 187 ;;; Lognormal distribution 188 189 (define (*makerandomlognormals mu sigma randoms) 190 (let ( 191 (normals (*makerandomnormals 0.0 1.0 randoms) ) 192 (nmu (log (* mu (/ mu (sqrt (+ (* sigma sigma) (* mu mu)))))) ) 193 (nsigma (sqrt (log (+ 1.0 (* sigma (/ sigma mu mu))))) ) ) 194 ; 195 (lambda () 196 (exp (+ nmu (* (normals) nsigma))))) ) 197 198 (define (makerandomlognormals #!key (mu 1.0) (sigma 1.0) (randoms (randomreal/current))) 199 (checknonzeroreal 'makerandomlognormals mu 'mu) 200 (checknonnegativereal 'makerandomlognormals sigma 'sigma) 201 (checkprocedure 'makerandomlognormals randoms 'randoms) 202 (values 203 (*makerandomlognormals mu sigma randoms) 204 (lambda () (values mu sigma randoms))) ) 205 206 ;;; Cauchy distribution 207 208 (define (*makerandomcauchys median sigma randoms) 209 (lambda () 210 (+ median (* sigma (tan (* *pi* ( (randoms) 0.5)))))) ) 211 212 (define (makerandomcauchys #!key (median 0.0) (sigma 1.0) (randoms (randomreal/current))) 213 (checkreal 'makerandomcauchys median 'median) 214 (checkpositivereal 'makerandomcauchys sigma 'sigma) 215 (checkprocedure 'makerandomcauchys randoms 'randoms) 216 (values 217 (*makerandomcauchys median sigma randoms) 218 (lambda () (values median sigma randoms))) ) 219 220 ;;; Gamma distribution 221 222 ;; "A Simple Method for Generating Gamma Variables", George Marsaglia & Wai Wan Tsang, 223 ;; ACM Transactions on Mathematical Software, Vol. 26, No. 3, September 2000, Pages 363 372. 224 225 (define (*makerandomgammas alpha theta randoms) 226 (if (= 1.0 alpha) 227 ;then special case 228 (if (= 1.0 theta) 229 (lambda () ( (log (randoms)))) 230 (lambda () (* theta ( (log (randoms))))) ) 231 ;else general case 232 (let ( 233 (normals (*makerandomnormals 0.0 1.0 randoms) ) 234 (uniforms 235 (if (< alpha 1.0) 236 (let ((alphainv (*reciprocal alpha))) 237 (lambda () (expt (randoms) alphainv) ) ) 238 randoms) ) ) 239 ; 240 (let* ( 241 (d ( (if (< alpha 1.0) (+ 1.0 alpha) alpha) *onethird*) ) 242 (c (*reciprocal (sqrt (* 9.0 d))) ) ) 243 ; 244 (lambda () 245 (* 246 theta 247 (let loop () 248 (let* ( 249 (x (normals) ) 250 (v (+ 1.0 (* c x)) ) ) 251 ; 252 (define (gamma?) 253 (let ( 254 (u (uniforms) ) 255 (x^2 (* x x) ) 256 (v^3 (* v v v) ) ) 257 ; 258 (or 259 (< u ( 1.0 (* 0.0331 x^2 x^2))) 260 (< (log u) (+ (* 0.5 x^2) (* d ( 1.0 (+ v^3 (log v^3)))))) ) ) ) 261 ; 262 (if (and (< 0.0 v) (gamma?)) 263 (* d v) 264 (loop) ) ) ) ) ) ) ) ) ) 265 266 (define (makerandomgammas #!key (alpha 1.0) (theta 1.0) (randoms (randomreal/current))) 267 (checkpositivereal 'makerandomgammas alpha 'alpha) 268 (checkpositivereal 'makerandomgammas theta 'theta) 269 (checkprocedure 'makerandomgammas randoms 'randoms) 270 (values 271 (*makerandomgammas alpha theta randoms) 272 (lambda () (values alpha theta randoms))) ) 273 274 ;;; Erlang distribution 275 276 (define (*makerandomerlangs alpha theta randoms) 277 (*makerandomgammas (exact>inexact alpha) (exact>inexact theta) randoms) ) 278 279 (define (makerandomerlangs #!key (alpha 1) (theta 1.0) (randoms (randomreal/current))) 280 (checkpositivereal 'makerandomerlangs alpha 'alpha) 281 (checkpositivereal 'makerandomerlangs theta 'theta) 282 (checkprocedure 'makerandomerlangs randoms 'randoms) 283 (values 284 (*makerandomerlangs alpha theta randoms) 285 (lambda () (values alpha theta randoms))) ) 286 287 ;;; Pareto distribution 288 289 (define (*makerandomparetos alpha xmin randoms) 290 (let ((gammas (*makerandomgammas alpha (*reciprocal xmin) randoms))) 291 (*makerandomexponentials 1.0 (lambda () (*reciprocal (+ xmin (gammas)))))) ) 292 293 (define (makerandomparetos #!key (alpha 1.0) (xmin 1.0) (randoms (randomreal/current))) 294 (checkpositivereal 'makerandomparetos alpha 'alpha) 295 (checkpositivereal 'makerandomparetos xmin 'xmin) 296 (checkprocedure 'makerandomparetos randoms 'randoms) 297 (values 298 (*makerandomparetos alpha xmin randoms) 299 (lambda () (values alpha xmin randoms))) ) 300 301 ;;; Levy distribution 302 303 ;; See Stable Distributions  John P. Nolan, Formula 1.12 304 305 (define (*makerandomlevys gamma delta randoms) 306 (if (and (= 1.0 gamma) (= 0.0 delta)) 307 (lambda () (let ((r (randoms))) (*reciprocal (* r r)))) 308 (lambda () (let ((r (randoms))) (+ delta (* gamma (*reciprocal (* r r))))))) ) 309 310 (define (makerandomlevys #!key (gamma 1.0) (delta 0.0) (randoms (randomreal/current))) 311 (checknonnegativereal 'makerandomlevys delta 'delta) 312 (checkpositivereal 'makerandomlevys gamma 'gamma) 313 (checkprocedure 'makerandomlevys randoms 'randoms) 314 (values 315 (*makerandomlevys gamma delta randoms) 316 (lambda () (values gamma delta randoms))) ) 317 318 ;;; Weibull distribution 319 320 (define (*makerandomweibulls shape scale randoms) 321 (let ((invscale (*reciprocal scale)) 322 (invshape (*reciprocal shape)) ) 323 (lambda () (expt (* invscale (log ( 1.0 (randoms)))) invshape)) ) ) 324 325 (define (makerandomweibulls #!key (shape 1.0) (scale 1.0) (randoms (randomreal/current))) 326 (checkpositivereal 'makerandomweibulls shape 'shape) 327 (checkpositivereal 'makerandomweibulls scale 'scale) 328 (checkprocedure 'makerandomweibulls randoms 'randoms) 329 (values 330 (*makerandomweibulls shape scale randoms) 331 (lambda () (values shape scale randoms))) ) 24 srfi27normals 25 srfi27exponentials 26 srfi27triangles 27 srfi27poissons 28 srfi27bernoullis 29 srfi27binomials 30 srfi27geometrics 31 srfi27lognormals 32 srfi27cauchys 33 srfi27gammas 34 srfi27erlangs 35 srfi27paretos 36 srfi27levys 37 srfi27weibulls) 332 38 333 39 ) ;module srfi27distributions 
release/4/srfi27/trunk/srfi27.meta
r34867 r34968 15 15 (numbers "2.8") 16 16 (mathh "3.2.0") 17 (dssslutils "2.1.0")17 #;(dssslutils "2.1.0") 18 18 #;(randombsd "0.2")) 19 19 (testdepends test) 20 20 (files 21 "srfi27.setup" "srfi27.meta" #;"srfi27.releaseinfo" 22 "srfi27.scm" 23 "srfi27implementation" 24 "randomsource.scm" 21 25 ;"bsdrnd.scm" 22 "compositerandomsource.scm" "compositeentropysource.scm" 23 "mrg32k3a.scm" "srfi27.meta" "srfi27uniformrandom.scm" "mwc.scm" 24 "srfi27.setup" "entropyprocedure.scm" "srfi27.releaseinfo" 25 "srfi27implementation" "srfi27numbers.scm" 26 "entropyunix.scm" "entropylinux.scm" "entropysource.scm" 27 "entropyclock.scm" "randomsource.scm" "srfi27vectorsupport.scm" 28 "srfi27distributions.scm" "srfi27vector.scm" "srfi27.scm" 29 "entropywindows.scm" "moa.scm" "entropysupport.scm" 26 "compositerandomsource.scm" 27 "mrg32k3a.scm" "mwc.scm" "moa.scm" 28 "srfi27numbers.scm" 29 "entropysource.scm" 30 "compositeentropysource.scm" 31 "entropyunix.scm" "entropylinux.scm" "entropywindows.scm" 32 "entropyclock.scm" "entropyprocedure.scm" "entropyport.scm" 33 "entropysupport.scm" 34 "srfi27uniformrandom.scm" 35 "srfi27distributions.scm" 36 "srfi27bernoullis.scm" 37 "srfi27binomials.scm" 38 "srfi27cauchys.scm" 39 "srfi27erlangs.scm" 40 "srfi27exponentials.scm" 41 "srfi27gammas.scm" 42 "srfi27geometrics.scm" 43 "srfi27levys.scm" 44 "srfi27lognormals.scm" 45 "srfi27normals.scm" 46 "srfi27paretos.scm" 47 "srfi27poissons.scm" 48 "srfi27triangles.scm" 49 "srfi27weibulls.scm" 50 "srfi27vector.scm" "srfi27vectorsupport.scm" 30 51 "sourceregistration.scm" "tests/testdiehard.scm" "tests/testconfidence.scm" 31 "tests/testmrg32k3a.scm" "tests/run.scm" "entropyport.scm") ) 52 "tests/testmrg32k3a.scm" 53 "tests/srfi27test.scm" "tests/run.scm") ) 
release/4/srfi27/trunk/srfi27.setup
r34967 r34968 124 124 #:compileoptions `(scrutinize ,@publoptn) ) 125 125 126 (setupsharedextensionmodule 'srfi27bernoullis (extensionversion "3.3.0") 127 #:inline? #t 128 #:types? #t 129 #:compileoptions `(scrutinize ,@publoptn) ) 130 131 (setupsharedextensionmodule 'srfi27binomials (extensionversion "3.3.0") 132 #:inline? #t 133 #:types? #t 134 #:compileoptions `(scrutinize ,@publoptn) ) 135 136 (setupsharedextensionmodule 'srfi27cauchys (extensionversion "3.3.0") 137 #:inline? #t 138 #:types? #t 139 #:compileoptions `(scrutinize ,@publoptn) ) 140 141 (setupsharedextensionmodule 'srfi27normals (extensionversion "3.3.0") 142 #:inline? #t 143 #:types? #t 144 #:compileoptions `(scrutinize ,@publoptn) ) 145 146 ;needs normals 147 (setupsharedextensionmodule 'srfi27gammas (extensionversion "3.3.0") 148 #:inline? #t 149 #:types? #t 150 #:compileoptions `(scrutinize ,@publoptn) ) 151 152 ;needs gammas 153 (setupsharedextensionmodule 'srfi27erlangs (extensionversion "3.3.0") 154 #:inline? #t 155 #:types? #t 156 #:compileoptions `(scrutinize ,@publoptn) ) 157 158 (setupsharedextensionmodule 'srfi27exponentials (extensionversion "3.3.0") 159 #:inline? #t 160 #:types? #t 161 #:compileoptions `(scrutinize ,@publoptn) ) 162 163 (setupsharedextensionmodule 'srfi27geometrics (extensionversion "3.3.0") 164 #:inline? #t 165 #:types? #t 166 #:compileoptions `(scrutinize ,@publoptn) ) 167 168 (setupsharedextensionmodule 'srfi27levys (extensionversion "3.3.0") 169 #:inline? #t 170 #:types? #t 171 #:compileoptions `(scrutinize ,@publoptn) ) 172 173 (setupsharedextensionmodule 'srfi27lognormals (extensionversion "3.3.0") 174 #:inline? #t 175 #:types? #t 176 #:compileoptions `(scrutinize ,@publoptn) ) 177 178 ;needs gammas exponentials 179 (setupsharedextensionmodule 'srfi27paretos (extensionversion "3.3.0") 180 #:inline? #t 181 #:types? #t 182 #:compileoptions `(scrutinize ,@publoptn) ) 183 184 (setupsharedextensionmodule 'srfi27poissons (extensionversion "3.3.0") 185 #:inline? #t 186 #:types? #t 187 #:compileoptions `(scrutinize ,@publoptn) ) 188 189 (setupsharedextensionmodule 'srfi27triangles (extensionversion "3.3.0") 190 #:inline? #t 191 #:types? #t 192 #:compileoptions `(scrutinize ,@publoptn) ) 193 194 (setupsharedextensionmodule 'srfi27weibulls (extensionversion "3.3.0") 195 #:inline? #t 196 #:types? #t 197 #:compileoptions `(scrutinize ,@publoptn) ) 198 126 199 (setupsharedextensionmodule 'srfi27distributions (extensionversion "3.3.0") 127 200 #:inline? #t
Note: See TracChangeset
for help on using the changeset viewer.