diff --git a/linear-least-squares.lisp b/linear-least-squares.lisp index 2eca08ce1adab8786564197491becaa414a00d05..ae3867b8cd4bdba32d8d16fc63c694b3e030f4f4 100644 --- a/linear-least-squares.lisp +++ b/linear-least-squares.lisp @@ -1,6 +1,6 @@ ;; Linear least squares, or linear regression ;; Liam Healy <2008-01-21 12:41:46EST linear-least-squares.lisp> -;; Time-stamp: <2008-02-02 21:01:35EST linear-least-squares.lisp> +;; Time-stamp: <2008-02-03 11:31:47EST linear-least-squares.lisp> ;; $Id: $ (in-package :gsl) @@ -290,13 +290,13 @@ (defun mv-linear-least-squares-data () "Generate data for second example in Section 36.5 of the GSL manual." - (rng-set *rng-mt19937* 0) - (loop for x from 1/10 below 2 by 1/10 - for xd = (coerce x 'double-float) - for y0 = (exp xd) - for sigma = (* 0.1d0 y0) - collect - (list xd (+ y0 (gaussian *rng-mt19937* sigma)) sigma))) + (letm ((rng (random-number-generator *mt19937* 0))) + (loop for x from 1/10 below 2 by 1/10 + for xd = (coerce x 'double-float) + for y0 = (exp xd) + for sigma = (* 0.1d0 y0) + collect + (list xd (+ y0 (gaussian rng sigma)) sigma)))) (defun mv-linear-least-squares-example (data) "Second example in Section 36.5 of the GSL manual." diff --git a/monte-carlo.lisp b/monte-carlo.lisp index 758304e3b1a30c7e74685f17048c3b10ba1fbd01..768ce2049719839853ff5515d9763fff25bf18d5 100644 --- a/monte-carlo.lisp +++ b/monte-carlo.lisp @@ -1,6 +1,6 @@ ;; Monte Carlo Integration ;; Liam Healy Sat Feb 3 2007 - 17:42 -;; Time-stamp: <2008-02-02 21:11:35EST monte-carlo.lisp> +;; Time-stamp: <2008-02-03 11:35:59EST monte-carlo.lisp> ;; $Id: $ (in-package :gsl) @@ -288,46 +288,25 @@ (def-mc-function monte-carlo-g 3) (defun random-walk-plain-example (&optional (nsamples 500000)) - (letm ((ws (monte-carlo-plain 3))) - (with-data (lower vector-double 3) - (with-data (upper vector-double 3) - (setf (data lower) #(0.0d0 0.0d0 0.0d0) - (data upper) (vector pi pi pi)) - (rng-set *rng-mt19937* 0) - (monte-carlo-integrate-plain - monte-carlo-g - lower upper - nsamples - *rng-mt19937* - ws))))) + (letm ((ws (monte-carlo-plain 3)) + (lower (vector-double #(0.0d0 0.0d0 0.0d0))) + (upper (vector-double (vector pi pi pi))) + (rng (random-number-generator *mt19937* 0))) + (monte-carlo-integrate-plain monte-carlo-g lower upper nsamples rng ws))))) (defun random-walk-miser-example (&optional (nsamples 500000)) - (letm ((ws (monte-carlo-miser 3))) - (with-data (lower vector-double 3) - (with-data (upper vector-double 3) - (setf (data lower) #(0.0d0 0.0d0 0.0d0) - (data upper) (vector pi pi pi)) - (rng-set *rng-mt19937* 0) - (monte-carlo-integrate-miser - monte-carlo-g - lower upper - nsamples - *rng-mt19937* - ws))))) + (letm ((ws (monte-carlo-miser 3)) + (lower (vector-double #(0.0d0 0.0d0 0.0d0))) + (upper (vector-double (vector pi pi pi))) + (rng (random-number-generator *mt19937* 0))) + (monte-carlo-integrate-miser monte-carlo-g lower upper nsamples rng ws))))) (defun random-walk-vegas-example (&optional (nsamples 500000)) - (letm ((ws (monte-carlo-vegas 3))) - (with-data (lower vector-double 3) - (with-data (upper vector-double 3) - (setf (data lower) #(0.0d0 0.0d0 0.0d0) - (data upper) (vector pi pi pi)) - (rng-set *rng-mt19937* 0) - (monte-carlo-integrate-vegas - monte-carlo-g - lower upper - nsamples - *rng-mt19937* - ws))))) + (letm ((ws (monte-carlo-vegas 3)) + (lower (vector-double #(0.0d0 0.0d0 0.0d0))) + (upper (vector-double (vector pi pi pi))) + (rng (random-number-generator *mt19937* 0))) + (monte-carlo-integrate-vegas monte-carlo-g lower upper nsamples rng ws))))) (lisp-unit:define-test monte-carlo (lisp-unit:assert-first-fp-equal diff --git a/random/bernoulli.lisp b/random/bernoulli.lisp index 69a80d792add25c8933220afedf4913acffbf780..ed9e8a0f4b1a867af72ffeea87304195d59e1044 100644 --- a/random/bernoulli.lisp +++ b/random/bernoulli.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: bernoulli.lisp -; description: Bernoulli distribution -; date: Sat Nov 25 2006 - 16:59 -; author: Liam M. Healy -; modified: Sat Nov 25 2006 - 17:01 -;******************************************************** -;;; $Id: $ +;; Bernoulli distribution +;; Liam Healy, Sat Nov 25 2006 - 16:59 +;; Time-stamp: <2008-02-03 11:16:35EST bernoulli.lisp> +;; $Id: $ (in-package :gsl) @@ -13,9 +9,9 @@ "gsl_ran_bernoulli" (((generator generator) :pointer) (p :double)) :c-return :uint - :documentation + :documentation ; FDL "Returns either 0 or 1, the result of a Bernoulli trial - with probability @var{p}. The probability distribution for + with probability p. The probability distribution for a Bernoulli trial is p(0) = 1 - p p(1) = p.") @@ -23,19 +19,19 @@ (defun-gsl bernoulli-pdf (k p) "gsl_ran_bernoulli_pdf" ((k :uint) (p :double)) :c-return :double - :documentation "The probability @math{p(k)} of obtaining - @var{k} from a Bernoulli distribution with probability parameter - @var{p}, using the formula given in #'bernoulli.") + :documentation ; FDL + "The probability p(k) of obtaining + k from a Bernoulli distribution with probability parameter + p, using the formula given in #'bernoulli.") ;;; Examples and unit test (lisp-unit:define-test bernoulli (lisp-unit:assert-equal '(0 1 1 0 1 1 0 0 0 0 0) - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 - collect - (bernoulli *rng-mt19937* 0.5d0)))) + collect + (bernoulli rng 0.5d0)))) (lisp-unit:assert-first-fp-equal "0.500000000000d+00" (bernoulli-pdf 0 0.5d0))) diff --git a/random/beta.lisp b/random/beta.lisp index 2e506106a3ce00f3213a5a4ccd221a1d507e5a04..d71a1747d5bd61b57a7b13a8d921a440be7b6f83 100644 --- a/random/beta.lisp +++ b/random/beta.lisp @@ -1,6 +1,6 @@ ;; Beta distribution ;; Liam Healy, Sat Sep 30 2006 -;; Time-stamp: <2008-02-02 23:04:13EST beta.lisp> +;; Time-stamp: <2008-02-03 09:52:01EST beta.lisp> ;; $Id: $ (in-package :gsl) @@ -10,42 +10,46 @@ "gsl_ran_beta" (((generator generator) :pointer) (a :double) (b :double)) :c-return :double - :documentation + :documentation ; FDL "A random variate from the beta distribution. The distribution function is p(x) dx = {\Gamma(a+b) \over \Gamma(a) \Gamma(b)} x^{a-1} (1-x)^{b-1} dx - @math{0 <= x <= 1}.") + 0 <= x <= 1.") (defun-gsl beta-pdf (x a b) "gsl_ran_beta_pdf" ((x :double) (a :double) (b :double)) :c-return :double - :documentation - "The probability density @math{p(x)} at @var{x} - for a beta distribution with parameters @var{a} and @var{b}, using the + :documentation ; FDL + "The probability density p(x) at x + for a beta distribution with parameters a and b, using the formula given in #'beta.") (defun-gsl beta-P (x a b) "gsl_cdf_beta_P" ((x :double) (a :double) (b :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{P(x)} for the beta distribution with parameters @var{a} and @var{b}.") + :documentation ; FDL + "The cumulative distribution functions + P(x) for the beta distribution with parameters a and b.") (defun-gsl beta-Q (x a b) "gsl_cdf_beta_Q" ((x :double) (a :double) (b :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{Q(x)} for the beta distribution with parameters @var{a} and @var{b}.") + :documentation ; FDL + "The cumulative distribution functions + Q(x) for the beta distribution with parameters a and b.") (defun-gsl beta-Pinv (P a b) "gsl_cdf_beta_Pinv" ((P :double) (a :double) (b :double)) :c-return :double - :documentation "The inverse cumulative distribution functions - @math{P(x)} for the beta distribution with parameters @var{a} and @var{b}.") + :documentation ; FDL + "The inverse cumulative distribution functions + P(x) for the beta distribution with parameters a and b.") (defun-gsl beta-Qinv (Q a b) "gsl_cdf_beta_Qinv" ((Q :double) (a :double) (b :double)) :c-return :double - :documentation "The inverse cumulative distribution functions - @math{Q(x)} for the beta distribution with parameters @var{a} and @var{b}.") + :documentation ; FDL + "The inverse cumulative distribution functions + Q(x) for the beta distribution with parameters a and b.") ;;; Examples and unit test (lisp-unit:define-test beta diff --git a/random/binomial.lisp b/random/binomial.lisp index d0a031778bde1006605c6727de17f88975960f11..f214d0266990e73af86e39ba4d7f40aa356b7ac7 100644 --- a/random/binomial.lisp +++ b/random/binomial.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: binomial.lisp -; description: Binomial distribution -; date: Sat Nov 25 2006 - 16:00 -; author: Liam M. Healy -; modified: Sat Nov 25 2006 - 18:08 -;******************************************************** -;;; $Id: $ +;; Binomial distribution +;; Liam Healy, Sat Nov 25 2006 - 16:00 +;; Time-stamp: <2008-02-03 11:25:08EST binomial.lisp> +;; $Id: $ (in-package :gsl) @@ -13,44 +9,44 @@ "gsl_ran_binomial" (((generator generator) :pointer) (p :double) (n :uint)) :c-return :uint - :documentation + :documentation ; FDL "A random integer from the binomial distribution, - the number of successes in @var{n} independent trials with probability - @var{p}. The probability distribution for binomial variates is, + the number of successes in n independent trials with probability + p. The probability distribution for binomial variates is, p(k) = {n! \over k! (n-k)!} p^k (1-p)^{n-k} - @math{0 <= k <= n}.") + 0 <= k <= n.") (defun-gsl binomial-pdf (k p n) "gsl_ran_binomial_pdf" ((k :uint) (p :double) (n :uint)) :c-return :double - :documentation - "The probability @math{p(k)} of obtaining @var{k} - from a binomial distribution with parameters @var{p} and @var{n}, using + :documentation ; FDL + "The probability p(k) of obtaining k + from a binomial distribution with parameters p and n, using the formula given in #'binomial.") (defun-gsl binomial-P (k p n) "gsl_cdf_binomial_P" ((k :uint) (p :double) (n :uint)) :c-return :double - :documentation "The cumulative distribution functions - @math{P(k)} for the Binomial distribution - with parameters @var{p} and @var{n}.") + :documentation ; FDL + "The cumulative distribution functions + P(k) for the Binomial distribution with parameters p and n.") (defun-gsl binomial-Q (k p n) "gsl_cdf_binomial_Q" ((k :uint) (p :double) (n :uint)) :c-return :double - :documentation "The cumulative distribution functions - @math{Q(k)} for the Binomial distribution - with parameters @var{p} and @var{n}.") + :documentation ; FDL + "The cumulative distribution functions + Q(k) for the Binomial distribution + with parameters p and n.") ;;; Examples and unit test (lisp-unit:define-test binomial (lisp-unit:assert-equal '(11 3 4 8 4 5 8 6 5 6 6) - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 - collect - (binomial *rng-mt19937* 0.4d0 12)))) + collect + (binomial rng 0.4d0 12)))) (lisp-unit:assert-first-fp-equal "0.227030335488d+00" (binomial-pdf 5 0.4d0 12)) diff --git a/random/cauchy.lisp b/random/cauchy.lisp index bb64c1bfbc07485bf245f894a39df146cbd40370..e05cdfbc3a143d6d06df562ed17383b62ac62125 100644 --- a/random/cauchy.lisp +++ b/random/cauchy.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: cauchy.lisp -; description: Cauchy distribution -; date: Sat Sep 30 2006 -; author: Liam M. Healy -; modified: Sat Sep 30 2006 - 19:15 -;******************************************************** -;;; $Id: $ +;; Cauchy distribution +;; Liam Healy, Sat Sep 30 2006 +;; Time-stamp: <2008-02-03 11:02:19EST cauchy.lisp> +;; $Id: $ (in-package :gsl) @@ -13,45 +9,49 @@ "gsl_ran_cauchy" (((generator generator) :pointer) (a :double)) :c-return :double - :documentation + :documentation ; FDL "A random variate from the Cauchy distribution with - scale parameter @var{a}. The probability distribution for Cauchy + scale parameter a. The probability distribution for Cauchy random variates is, p(x) dx = {1 \over a\pi (1 + (x/a)^2) } dx - for @math{x} in the range @math{-\infty} to @math{+\infty}. The Cauchy + for x in the range -\infty to +\infty. The Cauchy distribution is also known as the Lorentz distribution.") (defun-gsl cauchy-pdf (x a) "gsl_ran_cauchy_pdf" ((x :double) (a :double)) :c-return :double - :documentation - "The probability density @math{p(x)} at @var{x} - for a Cauchy distribution with scale parameter @var{a}, using the formula + :documentation ; FDL + "The probability density p(x) at x + for a Cauchy distribution with scale parameter a, using the formula given for #'cauchy.") (defun-gsl cauchy-P (x a) "gsl_cdf_cauchy_P" ((x :double) (a :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{P(x)} for the Cauchy distribution with scale parameter @var{a}.") + :documentation ; FDL + "The cumulative distribution functions + P(x) for the Cauchy distribution with scale parameter a.") (defun-gsl cauchy-Q (x a) "gsl_cdf_cauchy_Q" ((x :double) (a :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{Q(x)} for the Cauchy distribution with scale parameter @var{a}.") + :documentation ; FDL + "The cumulative distribution functions + Q(x) for the Cauchy distribution with scale parameter a.") (defun-gsl cauchy-Pinv (P a) "gsl_cdf_cauchy_Pinv" ((P :double) (a :double)) :c-return :double - :documentation "The inverse cumulative distribution functions - @math{P(x)} for the Cauchy distribution with scale parameter @var{a}.") + :documentation ; FDL + "The inverse cumulative distribution functions + P(x) for the Cauchy distribution with scale parameter a.") (defun-gsl cauchy-Qinv (Q a) "gsl_cdf_cauchy_Qinv" ((Q :double) (a :double)) :c-return :double - :documentation "The inverse cumulative distribution functions - @math{Q(x)} for the Cauchy distribution with scale parameter @var{a}.") + :documentation ; FDL + "The inverse cumulative distribution functions + Q(x) for the Cauchy distribution with scale parameter a.") ;;; Examples and unit test (lisp-unit:define-test cauchy @@ -61,11 +61,10 @@ "-0.134390491844d+01" "-0.103643632829d+02" "-0.790709314248d+02" "-0.106520710880d+02" "-0.939394824349d+01") (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 collect - (cauchy *rng-mt19937* 10.0d0))))) + (cauchy rng 10.0d0))))) (lisp-unit:assert-first-fp-equal "0.318309886184d-01" (cauchy-pdf 0.0d0 10.0d0)) diff --git a/random/chi-squared.lisp b/random/chi-squared.lisp index 7bc51cb976d4d3ced0beecb69a0de69b508bd119..1ae331d4a973f3711d3c4201374ffd1f23157a3f 100644 --- a/random/chi-squared.lisp +++ b/random/chi-squared.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: chi-squared.lisp -; description: Chi-squared distribution -; date: Sat Oct 7 2006 - 16:13 -; author: Liam M. Healy -; modified: Sun Oct 8 2006 - 16:41 -;******************************************************** -;;; $Id: $ +;; Chi-squared distribution +;; Liam Healy, Sat Oct 7 2006 - 16:13 +;; Time-stamp: <2008-02-03 10:23:21EST chi-squared.lisp> +;; $Id: $ (in-package :gsl) @@ -13,43 +9,47 @@ "gsl_ran_chisq" (((generator generator) :pointer) (nu :double)) :c-return :double - :documentation + :documentation ; FDL "A random variate from the chi-squared distribution - with @var{nu} degrees of freedom. The distribution function is + with nu degrees of freedom. The distribution function is p(x) dx = {1 \over 2 \Gamma(\nu/2) } (x/2)^{\nu/2 - 1} \exp(-x/2) dx - @math{x >= 0}. ") + x >= 0. ") (defun-gsl chi-squared-pdf (x nu) "gsl_ran_chisq_pdf" ((x :double) (nu :double)) :c-return :double - :documentation - "The probability density @math{p(x)} at @var{x} - for a chi-squared distribution with @var{nu} degrees of freedom, using + :documentation ; FDL + "The probability density p(x) at x + for a chi-squared distribution with nu degrees of freedom, using the formula given in #'chi-squared.") (defun-gsl chi-squared-P (x nu) "gsl_cdf_chisq_P" ((x :double) (nu :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{P(x)} for the chi-squared distribution with @var{nu} degrees of freedom.") + :documentation ; FDL + "The cumulative distribution functions + P(x) for the chi-squared distribution with nu degrees of freedom.") (defun-gsl chi-squared-Q (x nu) "gsl_cdf_chisq_Q" ((x :double) (nu :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{Q(x)} for the chi-squared distribution with @var{nu} degrees of freedom.") + :documentation ; FDL + "The cumulative distribution functions + Q(x) for the chi-squared distribution with nu degrees of freedom.") (defun-gsl chi-squared-Pinv (P nu) "gsl_cdf_chisq_Pinv" ((P :double) (nu :double)) :c-return :double - :documentation "The inverse cumulative distribution functions - @math{P(x)} for the chi-squared distribution with @var{nu} degrees of freedom.") + :documentation ; FDL + "The inverse cumulative distribution functions + P(x) for the chi-squared distribution with nu degrees of freedom.") (defun-gsl chi-squared-Qinv (Q nu) "gsl_cdf_chisq_Qinv" ((Q :double) (nu :double)) :c-return :double - :documentation "The inverse cumulative distribution functions - @math{Q(x)} for the chi-squared distribution with @var{nu} degrees of freedom.") + :documentation ; FDL + "The inverse cumulative distribution functions + Q(x) for the chi-squared distribution with nu degrees of freedom.") ;;; Examples and unit test (lisp-unit:define-test chi-squared @@ -59,11 +59,10 @@ "0.137566397216d+02" "0.113594580454d+02" "0.160682825906d+02" "0.129995142525d+02" "0.111473474244d+02") (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 collect - (chi-squared *rng-mt19937* 10.0d0))))) + (chi-squared rng 10.0d0))))) (lisp-unit:assert-first-fp-equal "0.439391289468d+00" (chi-squared-pdf 0.5d0 1.0d0)) diff --git a/random/dirichlet.lisp b/random/dirichlet.lisp index 618c46f2a496b99dfbd3100e2baedad70328909b..86bec3be5585d5bbc3602610da07706d3e230084 100644 --- a/random/dirichlet.lisp +++ b/random/dirichlet.lisp @@ -1,6 +1,6 @@ ;; Dirichlet distribution ;; Liam Healy, Sun Oct 29 2006 -;; Time-stamp: <2008-02-02 23:00:40EST dirichlet.lisp> +;; Time-stamp: <2008-02-03 09:49:17EST dirichlet.lisp> ;; $Id: $ (in-package :gsl) @@ -13,20 +13,20 @@ ;; theta had better be at least as long as alpha, or they'll be trouble ((gsl-array theta) :pointer)) :c-return :void - :documentation + :documentation ; FDL "An array of K=(length alpha) random variates from a Dirichlet distribution of order K-1. The distribution function is p(\theta_1,\ldots,\theta_K) \, d\theta_1 \cdots d\theta_K = {1 \over Z} \prod_{i=1}^{K} \theta_i^{\alpha_i - 1} \; \delta(1 -\sum_{i=1}^K \theta_i) d\theta_1 \cdots d\theta_K - @math{theta_i >= 0} and @math{alpha_i >= 0}. - The delta function ensures that @math{\sum \theta_i = 1}. - The normalization factor @math{Z} is + theta_i >= 0 and alpha_i >= 0. + The delta function ensures that \sum \theta_i = 1. + The normalization factor Z is Z = {\prod_{i=1}^K \Gamma(\alpha_i) \over \Gamma( \sum_{i=1}^K \alpha_i)} - The random variates are generated by sampling @var{K} values - from gamma distributions with parameters @math{a=alpha_i, b=1}, + The random variates are generated by sampling K values + from gamma distributions with parameters a=alpha_i, b=1, and renormalizing. - See A.M. Law, W.D. Kelton, @cite{Simulation Modeling and Analysis} + See A.M. Law, W.D. Kelton, \"Simulation Modeling and Analysis\" (1991).") (defun-gsl dirichlet-pdf (alpha theta) @@ -36,10 +36,10 @@ ;; theta had better be at least as long as alpha, or they'll be trouble ((gsl-array theta) :pointer)) :c-return :double - :documentation - "The probability density @math{p(\theta_1, ... , \theta_K)} - at @var{theta}[@var{K}] for a Dirichlet distribution with parameters - @var{alpha}[@var{K}], using the formula given for #'dirichlet.") + :documentation ; FDL + "The probability density p(\theta_1, ... , \theta_K) + at theta[K] for a Dirichlet distribution with parameters + alpha[K], using the formula given for #'dirichlet.") (defun-gsl dirichlet-log-pdf (alpha theta) "gsl_ran_dirichlet_lnpdf" @@ -48,11 +48,11 @@ ;; theta had better be at least as long as alpha, or they'll be trouble ((gsl-array theta) :pointer)) :c-return :double - :documentation + :documentation ; FDL "The logarithm of the probability density - @math{p(\theta_1, ... , \theta_K)} + p(\theta_1, ... , \theta_K) for a Dirichlet distribution with parameters - @var{alpha}[@var{K}].") + alpha[K].") ;;; Examples and unit test (lisp-unit:define-test dirichlet diff --git a/random/discrete.lisp b/random/discrete.lisp index 064a6ceabcc79bfcee2b9506b6d3a74a3d39330c..763ca6ba8d3442b2d9da15679c75596482246586 100644 --- a/random/discrete.lisp +++ b/random/discrete.lisp @@ -1,6 +1,6 @@ ;; Discrete random variables ;; Liam Healy, Sat Nov 11 2006 - 21:51 -;; Time-stamp: <2008-02-02 20:59:58EST discrete.lisp> +;; Time-stamp: <2008-02-03 10:58:16EST discrete.lisp> ;; $Id: $ (in-package :gsl) @@ -18,8 +18,8 @@ (((dim0 probabilities) :size) ((gsl-array probabilities) :pointer)) :c-return :pointer :export nil - :index (letm discrete) - :documentation + :index (letm discrete-random) + :documentation ; FDL "A pointer to a structure that contains the lookup table for the discrete random number generator. The array probabilities contains the probabilities of the discrete events; these array elements must all be @@ -31,22 +31,15 @@ "gsl_ran_discrete_free" ((table :pointer)) :c-return :void :export nil - :index (letm discrete) - :documentation + :index (letm discrete-random) + :documentation ; FDL "De-allocates the lookup table created by #'discrete-preprocess.") -(export 'with-discrete-table) -(defmacro with-discrete-table ((probabilities table) &body body) - `(let ((,table (discrete-preprocess ,probabilities))) - (unwind-protect - (progn ,@body) - (discrete-free ,table)))) - (defun-gsl discrete (generator table) "gsl_ran_discrete" (((generator generator) :pointer) (table :pointer)) :c-return :size - :documentation + :documentation ; FDL "Generate discrete random numbers after running #'discrete-preprocess; the argument 'table is the value returned by #'discrete-preprocess.") @@ -54,29 +47,28 @@ "gsl_ran_discrete_pdf" ((k :size) (table :pointer)) :c-return :double - :documentation - "The probability @math{P[k]} of observing the variable @var{k}. - Since @math{P[k]} is not stored as part of the lookup table, it must be - recomputed; this computation takes @math{O(K)}, so if @var{K} is large - and you care about the original array @math{P[k]} used to create the - lookup table, then you should just keep this original array @math{P[k]} + :documentation ; FDL + "The probability P[k] of observing the variable k. + Since P[k] is not stored as part of the lookup table, it must be + recomputed; this computation takes O(K), so if K is large + and you care about the original array P[k] used to create the + lookup table, then you should just keep this original array P[k] around.") ;;; Examples and unit test (lisp-unit:define-test discrete (lisp-unit:assert-equal '(1 0 1 1 0 1 1 2 1 2 2) - (with-data (probabilities vector-double 3) - (setf (data probabilities) #(0.25d0 0.5d0 0.25d0)) + ;; Must have two letms because the vector value is not set until + ;; the body, but the discrete-random needs that set value. + (letm ((probabilities (vector-double #(0.25d0 0.5d0 0.25d0)))) (letm ((table (discrete-random probabilities)) (rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 collect - (discrete *rng-mt19937* table)))) + (discrete rng table)))) (lisp-unit:assert-first-fp-equal "0.500000000000d+00" - (with-data (probabilities vector-double 3) - (setf (data probabilities) - #(0.25d0 0.5d0 0.25d0)) + (letm ((probabilities (vector-double #(0.25d0 0.5d0 0.25d0)))) (letm ((table (discrete-random probabilities))) (discrete-pdf 1 table)))))) diff --git a/random/exponential-power.lisp b/random/exponential-power.lisp index 86001d83b466a9a31929bb2ada3926962def448e..4d69f7083b1c510577ea1b76000a7b32155d01e0 100644 --- a/random/exponential-power.lisp +++ b/random/exponential-power.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: exponential-power.lisp -; description: Exponential power distribution -; date: Sat Sep 30 2006 -; author: Liam M. Healy -; modified: Sat Sep 30 2006 - 19:03 -;******************************************************** -;;; $Id: $ +;; Exponential power distribution +;; Liam Healy, Sat Sep 30 2006 +;; Time-stamp: <2008-02-03 11:12:27EST exponential-power.lisp> +;; $Id: $ (in-package :gsl) @@ -13,37 +9,38 @@ "gsl_ran_exppow" (((generator generator) :pointer) (a :double) (b :double)) :c-return :double - :documentation + :documentation ; FDL "A random variate from the exponential power distribution - with scale parameter @var{a} and exponent @var{b}. The distribution is + with scale parameter a and exponent b. The distribution is p(x) dx = {1 \over 2 a \Gamma(1+1/b)} \exp(-|x/a|^b) dx - for @math{x >= 0}. For @math{b = 1} this reduces to the Laplace - distribution. For @math{b = 2} it has the same form as a gaussian - distribution, but with @c{$a = \sqrt{2} \sigma$} - @math{a = \sqrt@{2@} \sigma}.") + for x >= 0. For b = 1 this reduces to the Laplace + distribution. For b = 2 it has the same form as a gaussian + distribution, but with a = \sqrt{2} \sigma.") (defun-gsl exponential-power-pdf (x a b) "gsl_ran_exppow_pdf" ((x :double) (a :double) (b :double)) :c-return :double - :documentation - "The probability density @math{p(x)} at @var{x} - for an exponential power distribution with scale parameter @var{a} - and exponent @var{b}, using the formula given for #'exponential-power.") + :documentation ; FDL + "The probability density p(x) at x + for an exponential power distribution with scale parameter a + and exponent b, using the formula given for #'exponential-power.") (defun-gsl exponential-power-P (x a b) "gsl_cdf_exppow_P" ((x :double) (a :double) (b :double)) :c-return :double - :documentation "The cumulative distribution function - @math{P(x)}, for the exponential power distribution with - parameters @var{a} and @var{b}.") + :documentation ; FDL + "The cumulative distribution function + P(x), for the exponential power distribution with + parameters a and b.") (defun-gsl exponential-power-Q (x a b) "gsl_cdf_exppow_Q" ((x :double) (a :double) (b :double)) :c-return :double - :documentation "The cumulative distribution functions @math{Q(x)} + :documentation ; FDL + "The cumulative distribution functions Q(x) for the exponential power distribution with - parameters @var{a} and @var{b}.") + parameters a and b.") ;;; Examples and unit test (lisp-unit:define-test exponential-power @@ -53,11 +50,10 @@ "-0.169473362899d+01" "-0.480323610806d+00" "-0.276417363499d-01" "0.631839185605d+00" "-0.124788752274d-01") (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 collect - (exponential-power *rng-mt19937* 1.0d0 2.0d0))))) + (exponential-power rng 1.0d0 2.0d0))))) (lisp-unit:assert-first-fp-equal "0.564189583548d+00" (exponential-power-pdf 0.0d0 1.0d0 2.0d0)) diff --git a/random/exponential.lisp b/random/exponential.lisp index 13b5b2480176ea9515a6fcadcb32eb29bba090ce..b8582069bc7a336dd06352c98716c78cf5d09b04 100644 --- a/random/exponential.lisp +++ b/random/exponential.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: exponential.lisp -; description: Exponential distribution -; date: Sat Sep 2 2006 - 19:04 -; author: Liam M. Healy -; modified: Mon Sep 11 2006 - 11:58 -;******************************************************** -;;; $Id: $ +;; Exponential distribution +;; Liam Healy, Sat Sep 2 2006 - 19:04 +;; Time-stamp: <2008-02-03 10:42:55EST exponential.lisp> +;; $Id: $ (in-package :gsl) @@ -13,43 +9,47 @@ "gsl_ran_exponential" (((generator generator) :pointer) (mu :double)) :c-return :double - :documentation + :documentation ; FDL "A random variate from the exponential distribution - with mean @var{mu}. The distribution is + with mean mu. The distribution is p(x) dx = {1 \over \mu} \exp(-x/\mu) dx - @math{x >= 0}.") + x >= 0.") (defun-gsl exponential-pdf (x mu) "gsl_ran_exponential_pdf" ((x :double) (mu :double)) :c-return :double - :documentation - "The probability density @math{p(x)} at @var{x} - for an exponential distribution with mean @var{mu}, using the formula + :documentation ; FDL + "The probability density p(x) at x + for an exponential distribution with mean mu, using the formula given for exponential.") (defun-gsl exponential-P (x mu) "gsl_cdf_exponential_P" ((x :double) (mu :double)) :c-return :double - :documentation "The cumulative distribution function - @math{P(x)} for the exponential distribution with mean @var{mu}.") + :documentation ; FDL + "The cumulative distribution function + P(x) for the exponential distribution with mean mu.") (defun-gsl exponential-Q (x mu) "gsl_cdf_exponential_Q" ((x :double) (mu :double)) :c-return :double - :documentation "The cumulative distribution function - @math{Q(x)} for the exponential distribution with mean @var{mu}.") + :documentation ; FDL + "The cumulative distribution function + Q(x) for the exponential distribution with mean mu.") (defun-gsl exponential-Pinv (P mu) "gsl_cdf_exponential_Pinv" ((P :double) (mu :double)) :c-return :double - :documentation "The inverse cumulative distribution function - @math{P(x)} for the exponential distribution with mean @var{mu}.") + :documentation ; FDL + "The inverse cumulative distribution function + P(x) for the exponential distribution with mean mu.") (defun-gsl exponential-Qinv (Q mu) "gsl_cdf_exponential_Qinv" ((Q :double) (mu :double)) :c-return :double - :documentation "The inverse cumulative distribution function - @math{Q(x)} for the exponential distribution with mean @var{mu}.") + :documentation ; FDL + "The inverse cumulative distribution function + Q(x) for the exponential distribution with mean mu.") ;;; Examples and unit test (lisp-unit:define-test exponential @@ -59,11 +59,10 @@ "0.434536244968d+00" "0.295303920905d+01" "0.616105293907d+01" "0.301168633354d+01" "0.274510798194d+01") (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 collect - (exponential *rng-mt19937* 10.0d0))))) + (exponential rng 10.0d0))))) (lisp-unit:assert-first-fp-equal "0.100000000000d+00" (exponential-pdf 0.0d0 10.0d0)) @@ -79,3 +78,4 @@ (lisp-unit:assert-first-fp-equal "0.100000000000d+01" (exponential-qinv 0.6065306597126334d0 2.0d0))) + diff --git a/random/fdist.lisp b/random/fdist.lisp index 979b6988731816db47a43f69b498973d74642f67..a15dee02879d734bf98b227f543a02259110c19c 100644 --- a/random/fdist.lisp +++ b/random/fdist.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: fdist.lisp -; description: Fdist distribution -; date: Sat Sep 30 2006 -; author: Liam M. Healy -; modified: Sun Oct 8 2006 - 16:31 -;******************************************************** -;;; $Id: $ +;; Fdist distribution +;; Liam Healy, Sat Sep 30 2006 +;; Time-stamp: <2008-02-03 10:17:28EST fdist.lisp> +;; $Id: $ (in-package :gsl) @@ -13,51 +9,55 @@ "gsl_ran_fdist" (((generator generator) :pointer) (nu1 :double) (nu2 :double)) :c-return :double - :documentation - "A random variate from the F-distribution with degrees of freedom @var{nu1} - and @var{nu2}. The distribution function is + :documentation ; FDL + "A random variate from the F-distribution with degrees of freedom nu1 + and nu2. The distribution function is p(x) dx = { \Gamma((\nu_1 + \nu_2)/2) \over \Gamma(\nu_1/2) \Gamma(\nu_2/2) } \nu_1^{\nu_1/2} \nu_2^{\nu_2/2} x^{\nu_1/2 - 1} (\nu_2 + \nu_1 x)^{-\nu_1/2 -\nu_2/2} - @math{x >= 0}.") + for x >= 0.") (defun-gsl fdist-pdf (x nu1 nu2) "gsl_ran_fdist_pdf" ((x :double) (nu1 :double) (nu2 :double)) :c-return :double - :documentation - "The probability density @math{p(x)} at @var{x} - for an F-distribution with @var{nu1} and @var{nu2} degrees of freedom, + :documentation ; FDL + "The probability density p(x) at x + for an F-distribution with nu1 and nu2 degrees of freedom, using the formula given #'fdist.") (defun-gsl fdist-P (x nu1 nu2) "gsl_cdf_fdist_P" ((x :double) (nu1 :double) (nu2 :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{P(x)} for the fdist distribution with - @var{nu1} and @var{nu2} degrees of freedom.") + :documentation ; FDL + "The cumulative distribution functions + P(x) for the fdist distribution with + nu1 and nu2 degrees of freedom.") (defun-gsl fdist-Q (x nu1 nu2) "gsl_cdf_fdist_Q" ((x :double) (nu1 :double) (nu2 :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{Q(x)} for the fdist distribution with - @var{nu1} and @var{nu2} degrees of freedom.") + :documentation ; FDL + "The cumulative distribution functions + Q(x) for the fdist distribution with + nu1 and nu2 degrees of freedom.") (defun-gsl fdist-Pinv (P nu1 nu2) "gsl_cdf_fdist_Pinv" ((P :double) (nu1 :double) (nu2 :double)) :c-return :double - :documentation "The inverse cumulative distribution functions - @math{P(x)} for the fdist distribution with - @var{nu1} and @var{nu2} degrees of freedom.") + :documentation ; FDL + "The inverse cumulative distribution functions + P(x) for the fdist distribution with + nu1 and nu2 degrees of freedom.") (defun-gsl fdist-Qinv (Q nu1 nu2) "gsl_cdf_fdist_Qinv" ((Q :double) (nu1 :double) (nu2 :double)) :c-return :double - :documentation "The inverse cumulative distribution functions - @math{Q(x)} for the fdist distribution with - @var{nu1} and @var{nu2} degrees of freedom.") + :documentation ; FDL + "The inverse cumulative distribution functions + Q(x) for the fdist distribution with + nu1 and nu2 degrees of freedom.") ;;; Examples and unit test (lisp-unit:define-test fdist @@ -67,11 +67,10 @@ "0.635729092565d-01" "0.477966365217d+00" "0.347211676079d-01" "0.486974823041d+00" "0.253179451696d+01") (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 collect - (fdist *rng-mt19937* 1.0d0 2.0d0))))) + (fdist rng 1.0d0 2.0d0))))) (lisp-unit:assert-first-fp-equal "0.159471988462d+00" (fdist-pdf 1.2d0 1.0d0 2.0d0)) diff --git a/random/flat.lisp b/random/flat.lisp index d19686ad87f4c9648464977d021885f65643a492..4cc16d5b294de1eeaa3f91b754d2417951ea9158 100644 --- a/random/flat.lisp +++ b/random/flat.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: flat.lisp -; description: Flat distribution -; date: Oct 7 2006 -; author: Liam M. Healy -; modified: Sat Oct 7 2006 - 16:13 -;******************************************************** -;;; $Id: $ +;; Flat distribution +;; Liam Healy, Oct 7 2006 +;; Time-stamp: <2008-02-03 11:18:30EST flat.lisp> +;; $Id: $ (in-package :gsl) @@ -13,43 +9,47 @@ "gsl_ran_flat" (((generator generator) :pointer) (a :double) (b :double)) :c-return :double - :documentation + :documentation ; FDL "A random variate from the flat (uniform) - distribution from @var{a} to @var{b}. The distribution is, + distribution from a to b. The distribution is p(x) dx = {1 \over (b-a)} dx - if @c{$a \le x < b$} @math{a <= x < b} and 0 otherwise.") + if a <= x < b, and 0 otherwise.") (defun-gsl flat-pdf (x a b) "gsl_ran_flat_pdf" ((x :double) (a :double) (b :double)) :c-return :double - :documentation - "The probability density @math{p(x)} at @var{x} - for a uniform distribution from @var{a} to @var{b}, using the formula + :documentation ; FDL + "The probability density p(x) at x + for a uniform distribution from a to b, using the formula given for #'flat.") (defun-gsl flat-P (x a b) "gsl_cdf_flat_P" ((x :double) (a :double) (b :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{P(x)} for a uniform distribution from @var{a} to @var{b}.") + :documentation ; FDL + "The cumulative distribution functions + P(x) for a uniform distribution from a to b.") (defun-gsl flat-Q (x a b) "gsl_cdf_flat_Q" ((x :double) (a :double) (b :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{Q(x)} for a uniform distribution from @var{a} to @var{b}.") + :documentation ; FDL + "The cumulative distribution functions + Q(x) for a uniform distribution from a to b.") (defun-gsl flat-Pinv (P a b) "gsl_cdf_flat_Pinv" ((P :double) (a :double) (b :double)) :c-return :double - :documentation "The inverse cumulative distribution functions - @math{P(x)} for a uniform distribution from @var{a} to @var{b}.") + :documentation ; FDL + "The inverse cumulative distribution functions + P(x) for a uniform distribution from a to b.") (defun-gsl flat-Qinv (Q a b) "gsl_cdf_flat_Qinv" ((Q :double) (a :double) (b :double)) :c-return :double - :documentation "The inverse cumulative distribution functions - @math{Q(x)} for a uniform distribution from @var{a} to @var{b}.") + :documentation ; FDL + "The inverse cumulative distribution functions + Q(x) for a uniform distribution from a to b.") ;;; Examples and unit test (lisp-unit:define-test flat @@ -59,11 +59,10 @@ "0.195747695654d+01" "0.174430534313d+01" "0.154004365834d+01" "0.173995298147d+01" "0.175994379818d+01") (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 collect - (flat *rng-mt19937* 1.0d0 2.0d0))))) + (flat rng 1.0d0 2.0d0))))) (lisp-unit:assert-first-fp-equal "0.100000000000d+01" (flat-pdf 1.2d0 1.0d0 2.0d0)) diff --git a/random/gamma.lisp b/random/gamma.lisp index a439c6a6b6ad2f670b4239e30032db9d9b7a1e8b..0df209fee850fc0e999d509fc0d608cb010cc773 100644 --- a/random/gamma.lisp +++ b/random/gamma.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: gamma.lisp -; description: Gamma distribution -; date: Sat Sep 30 2006 -; author: Liam M. Healy -; modified: Sat Oct 7 2006 - 11:53 -;******************************************************** -;;; $Id: $ +;; Gamma distribution +;; Liam Healy, Sat Sep 30 2006 +;; Time-stamp: <2008-02-03 10:43:54EST gamma.lisp> +;; $Id: $ (in-package :gsl) @@ -14,11 +10,11 @@ "gsl_ran_gamma" (((generator generator) :pointer) (a :double) (b :double)) :c-return :double - :documentation + :documentation ; FDL "A random variate from the gamma distribution. The distribution function is p(x) dx = {1 \over \Gamma(a) b^a} x^{a-1} e^{-x/b} dx - for @math{x > 0}. The gamma distribution with an integer parameter @var{a} + for x > 0. The gamma distribution with an integer parameter a is known as the Erlang distribution. The variates are computed using the algorithms from Knuth (vol 2).") @@ -26,40 +22,44 @@ "gsl_ran_gamma_mt" (((generator generator) :pointer) (a :double) (b :double)) :c-return :double - :documentation + :documentation ; FDL "A gamma variate using the Marsaglia-Tsang fast gamma method.") (defun-gsl gamma-pdf (x a b) "gsl_ran_gamma_pdf" ((x :double) (a :double) (b :double)) :c-return :double - :documentation - "The probability density @math{p(x)} at @var{x} - for a gamma distribution with parameters @var{a} and @var{b}, using the + :documentation ; FDL + "The probability density p(x) at x + for a gamma distribution with parameters a and b, using the formula given in #'gamma.") (defun-gsl gamma-P (x a b) "gsl_cdf_gamma_P" ((x :double) (a :double) (b :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{P(x)} for the Gamma distribution with parameters @var{a} and @var{b}.") + :documentation ; FDL + "The cumulative distribution functions + P(x) for the Gamma distribution with parameters a and b.") (defun-gsl gamma-Q (x a b) "gsl_cdf_gamma_Q" ((x :double) (a :double) (b :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{Q(x)} for the Gamma distribution with parameters @var{a} and @var{b}.") + :documentation ; FDL + "The cumulative distribution functions + Q(x) for the Gamma distribution with parameters a and b.") (defun-gsl gamma-Pinv (P a b) "gsl_cdf_gamma_Pinv" ((P :double) (a :double) (b :double)) :c-return :double - :documentation "The inverse cumulative distribution functions - @math{P(x)} for the Gamma distribution with parameters @var{a} and @var{b}.") + :documentation ; FDL + "The inverse cumulative distribution functions + P(x) for the Gamma distribution with parameters a and b.") (defun-gsl gamma-Qinv (Q a b) "gsl_cdf_gamma_Qinv" ((Q :double) (a :double) (b :double)) :c-return :double - :documentation "The inverse cumulative distribution functions - @math{Q(x)} for the Gamma distribution with parameters @var{a} and @var{b}.") + :documentation ; FDL + "The inverse cumulative distribution functions + Q(x) for the Gamma distribution with parameters a and b.") ;;; Examples and unit test (lisp-unit:define-test gamma-randist @@ -69,22 +69,20 @@ "0.869072489937d-01" "0.590607841809d+00" "0.123221058781d+01" "0.602337266708d+00" "0.549021596387d+00") (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 collect - (gamma-rd *rng-mt19937* 1.0d0 2.0d0))))) + (gamma-rd rng 1.0d0 2.0d0))))) (lisp-unit:assert-equal '("0.260013787613d+01" "0.252226669542d+01" "0.773142209213d+01" "0.422727992649d+01" "0.951930434749d-01" "0.571092010687d+00" "0.891063771946d+00" "0.826322120255d+00" "0.318306657206d+01" "0.380840036132d-02" "0.103201173341d+01") (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 collect - (gamma-mt *rng-mt19937* 1.0d0 2.0d0))))) + (gamma-mt rng 1.0d0 2.0d0))))) (lisp-unit:assert-first-fp-equal "0.475614712250d+00" (gamma-pdf 0.1d0 1.0d0 2.0d0)) diff --git a/random/gaussian-bivariate.lisp b/random/gaussian-bivariate.lisp index e4f1e1ff7cf76b968b2c88ea3d6ae0dc0803ef4f..432ccb5ed79a75ca780e5bf9a0cd31c23036bd2a 100644 --- a/random/gaussian-bivariate.lisp +++ b/random/gaussian-bivariate.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: gaussian-bivariate.lisp -; description: Gaussian bivariate distribution -; date: Sat Sep 2 2006 - 16:32 -; author: Liam M. Healy -; modified: Sat Sep 2 2006 - 21:34 -;******************************************************** -;;; $Id: $ +;; Gaussian bivariate distribution +;; Liam Healy, Sat Sep 2 2006 - 16:32 +;; Time-stamp: <2008-02-03 10:30:48EST gaussian-bivariate.lisp> +;; $Id: $ (in-package :gsl) @@ -14,27 +10,27 @@ (((generator generator) :pointer) (sigma-x :double) (sigma-y :double) (rho :double) (x :double) (y :double)) :c-return :void - :documentation + :documentation ; FDL "Generate a pair of correlated Gaussian variates, with - mean zero, correlation coefficient @var{rho} and standard deviations - @var{sigma_x} and @var{sigma_y} in the @math{x} and @math{y} directions. + mean zero, correlation coefficient rho and standard deviations + sigma_x and sigma_y in the x and y directions. The probability distribution for bivariate Gaussian random variates is, p(x,y) dx dy = {1 \over 2 \pi \sigma_x \sigma_y \sqrt{1-\rho^2}} \exp \left(-{(x^2/\sigma_x^2 + y^2/\sigma_y^2 - 2 \rho x y/(\sigma_x\sigma_y)) \over 2(1-\rho^2)}\right) dx dy - for @math{x,y} in the range @math{-\infty} to @math{+\infty}. The - correlation coefficient @var{rho} should lie between @math{1} and @math{-1}.") + for x,y in the range -\infty to +\infty. The + correlation coefficient rho should lie between 1 and -1.") (defun-gsl bivariate-gaussian-pdf (x y sigma-x sigma-y rho) "gsl_ran_bivariate_gaussian_pdf" ((x :double) (y :double) (sigma-x :double) (sigma-y :double) (rho :double)) :c-return :double - :documentation - "The probability density @math{p(x,y)} at - (@var{x},@var{y}) for a bivariate Gaussian distribution with standard - deviations @var{sigma_x}, @var{sigma_y} and correlation coefficient - @var{rho}, using the formula given for bivariate-gaussian.") + :documentation ; FDL + "The probability density p(x,y) at + (x,y) for a bivariate Gaussian distribution with standard + deviations sigma_x, sigma_y and correlation coefficient + rho, using the formula given for bivariate-gaussian.") ;;; Examples and unit test (lisp-unit:define-test gaussian-bivariate @@ -44,11 +40,10 @@ "-0.544622941217d+00" "-0.659202684161d+00" "-0.110295166108d+00" "0.179318404121d+00" "0.210251049803d+01") (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 collect - (bivariate-gaussian *rng-mt19937* 1.0d0 0.75d0 0.25d0))))) + (bivariate-gaussian rng 1.0d0 0.75d0 0.25d0))))) (lisp-unit:assert-first-fp-equal "0.554826555797d+00" (bivariate-gaussian-pdf 0.25d0 0.5d0 0.25d0 diff --git a/random/gaussian-tail.lisp b/random/gaussian-tail.lisp index 583d1ff69d7aad777d706b9e5919052b7239b84a..3e1715aea7d17be655c657c27a8047514ce0f378 100644 --- a/random/gaussian-tail.lisp +++ b/random/gaussian-tail.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: gaussian-tail.lisp -; description: Gaussian tail distribution -; date: Mon Aug 21 2006 - 21:52 -; author: Liam M. Healy -; modified: Mon Aug 21 2006 - 22:20 -;******************************************************** -;;; $Id: $ +;; Gaussian tail distribution +;; Liam Healy, Mon Aug 21 2006 - 21:52 +;; Time-stamp: <2008-02-03 10:15:01EST gaussian-tail.lisp> +;; $Id: $ (in-package :gsl) @@ -13,38 +9,38 @@ "gsl_ran_gaussian_tail" (((generator generator) :pointer) (a :double) (sigma :double)) :c-return :double - :documentation + :documentation ; FDL "Random variates from the upper tail of a Gaussian - distribution with standard deviation @var{sigma}. The values returned - are larger than the lower limit @var{a}, which must be positive. The + distribution with standard deviation sigma. The values returned + are larger than the lower limit a, which must be positive. The method is based on Marsaglia's famous rectangle-wedge-tail algorithm (Ann. Math. Stat. 32, 894--899 (1961)), with this aspect explained in Knuth, v2, 3rd ed, p139,586 (exercise 11). The probability distribution for Gaussian tail random variates is, p(x) dx = {1 \over N(a;\sigma) \sqrt{2 \pi \sigma^2}} \exp (- x^2 / 2\sigma^2) dx - for @math{x > a} where @math{N(a;\sigma)} is the normalization constant, + for x > a where N(a;\sigma) is the normalization constant, N(a;\sigma) = (1/2) erfc(a / sqrt(2 sigma^2)).") (defun-gsl gaussian-tail-pdf (x a sigma) "gsl_ran_gaussian_tail_pdf" ((x :double) (a :double) (sigma :double)) :c-return :double - :documentation - "The probability density @math{p(x)} at @var{x} - for a Gaussian tail distribution with standard deviation @var{sigma} and - lower limit @var{a}, using the formula given for gaussian-tail.") + :documentation ; FDL + "The probability density p(x) at x + for a Gaussian tail distribution with standard deviation sigma and + lower limit a, using the formula given for gaussian-tail.") (defun-gsl ugaussian-tail (generator a) "gsl_ran_ugaussian_tail" (((generator generator) :pointer) (a :double)) :c-return :double - :documentation + :documentation ; FDL "Equivalent to gaussian-tail with sigma=1.") (defun-gsl ugaussian-tail-pdf (x a) "gsl_ran_ugaussian_tail_pdf" ((x :double) (a :double)) :c-return :double - :documentation + :documentation ; FDL "Equivalent to gaussian-tail-pdf with sigma=1.") ;;; Examples and unit test @@ -55,11 +51,10 @@ "0.534422862873d+02" "0.518376171418d+02" "0.530010742143d+02" "0.521497741699d+02" "0.501157244350d+02") (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 collect - (gaussian-tail *rng-mt19937* 50.0d0 10.0d0))))) + (gaussian-tail rng 50.0d0 10.0d0))))) (lisp-unit:assert-first-fp-equal "0.187022708773d+00" (gaussian-tail-pdf 52.0d0 50.0d0 10.0d0)) @@ -69,11 +64,10 @@ "0.534422862873d+01" "0.518376171418d+01" "0.530010742143d+01" "0.521497741699d+01" "0.501157244350d+01") (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 collect - (ugaussian-tail *rng-mt19937* 5.0d0))))) + (ugaussian-tail rng 5.0d0))))) (lisp-unit:assert-first-fp-equal "0.187022708773d+01" (ugaussian-tail-pdf 5.2d0 5.0d0))) diff --git a/random/gaussian.lisp b/random/gaussian.lisp index c43bcbf2dec7c54b780140f09219f96fd1d138e9..6f0bc6e55e415902c5ca15ecb0a69d98556f3aff 100644 --- a/random/gaussian.lisp +++ b/random/gaussian.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: gaussian.lisp -; description: Gaussian distribution -; date: Sun Jul 16 2006 - 22:09 -; author: Liam M. Healy -; modified: Sun Jul 16 2006 - 22:42 -;******************************************************** -;;; $Id: $ +;; Gaussian distribution +;; Liam Healy, Sun Jul 16 2006 - 22:09 +;; Time-stamp: <2008-02-03 10:21:12EST gaussian.lisp> +;; $Id: $ (in-package :gsl) @@ -13,29 +9,29 @@ "gsl_ran_gaussian" (((generator generator) :pointer) (sigma :double)) :c-return :double - :documentation + :documentation ; FDL "A Gaussian random variate, with mean zero and - standard deviation @var{sigma}. The probability distribution for + standard deviation sigma. The probability distribution for Gaussian random variates is p(x) dx = {1 \over \sqrt{2 \pi \sigma^2}} \exp (-x^2 / 2\sigma^2) dx - for @math{x} in the range @math{-\infty} to @math{+\infty}. Use the - transformation @math{z = \mu + x} on the numbers returned by - @code{gsl_ran_gaussian} to obtain a Gaussian distribution with mean - @math{\mu}. This function uses the Box-Mueller algorithm which requires two - calls to the random number generator @var{r}.") + for x in the range -\infty to +\infty. Use the + transformation z = \mu + x on the numbers returned by + #'gaussian to obtain a Gaussian distribution with mean + mu. This function uses the Box-Mueller algorithm which requires two + calls to the random number generator r.") (defun-gsl gaussian-pdf (x sigma) "gsl_ran_gaussian_pdf" ((x :double) (sigma :double)) :c-return :double - :documentation - "Compute the probability density @math{p(x)} at @var{x} - for a Gaussian distribution with standard deviation @var{sigma}.") + :documentation ; FDL + "Compute the probability density p(x) at x + for a Gaussian distribution with standard deviation sigma.") (defun-gsl gaussian-ziggurat (generator sigma) "gsl_ran_gaussian_ziggurat" (((generator generator) :pointer) (sigma :double)) :c-return :double - :documentation + :documentation ; FDL "Compute a Gaussian random variate using the alternative Marsaglia-Tsang ziggurat method. The Ziggurat algorithm is the fastest available algorithm in most cases.") @@ -44,89 +40,89 @@ "gsl_ran_gaussian_ratio_method" (((generator generator) :pointer) (sigma :double)) :c-return :double - :documentation + :documentation ; FDL "Compute a Gaussian random variate using the Kinderman-Monahan-Leva ratio method.") (defun-gsl ugaussian (generator) "gsl_ran_ugaussian" (((generator generator) :pointer)) :c-return :double - :documentation + :documentation ; FDL "Compute results for the unit Gaussian distribution, equivalent to the #'gaussian with a standard deviation of one, - @var{sigma} = 1.") + sigma = 1.") (defun-gsl ugaussian-pdf (x) "gsl_ran_ugaussian_pdf" ((x :double)) :c-return :double - :documentation + :documentation ; FDL "Compute results for the unit Gaussian distribution, equivalent to the #'gaussian-pdf with a standard deviation of one, - @var{sigma} = 1.") + sigma = 1.") (defun-gsl ugaussian-ratio-method (generator) "gsl_ran_ugaussian_ratio_method" (((generator generator) :pointer)) :c-return :double - :documentation + :documentation ; FDL "Compute results for the unit Gaussian distribution, equivalent to the #'gaussian-ration-method with a - standard deviation of one, @var{sigma} = 1.") + standard deviation of one, sigma = 1.") (defun-gsl gaussian-P (x sigma) "gsl_cdf_gaussian_P" ((x :double) (sigma :double)) :c-return :double - :documentation - "The cumulative distribution function @math{P(x)} for the Gaussian - distribution with standard deviation @var{sigma}.") + :documentation ; FDL + "The cumulative distribution function P(x) for the Gaussian + distribution with standard deviation sigma.") (defun-gsl gaussian-Q (x sigma) "gsl_cdf_gaussian_Q" ((x :double) (sigma :double)) :c-return :double - :documentation - "The cumulative distribution function @math{Q(x)} for the Gaussian - distribution with standard deviation @var{sigma}.") + :documentation ; FDL + "The cumulative distribution function Q(x) for the Gaussian + distribution with standard deviation sigma.") (defun-gsl gaussian-Pinv (P sigma) "gsl_cdf_gaussian_Pinv" ((P :double) (sigma :double)) :c-return :double - :documentation - "The inverse cumulative distribution function @math{P(x)} for the Gaussian - distribution with standard deviation @var{sigma}.") + :documentation ; FDL + "The inverse cumulative distribution function P(x) for the Gaussian + distribution with standard deviation sigma.") (defun-gsl gaussian-Qinv (Q sigma) "gsl_cdf_gaussian_Qinv" ((Q :double) (sigma :double)) :c-return :double - :documentation - "The inverse cumulative distribution function @math{Q(x)} for the Gaussian - distribution with standard deviation @var{sigma}.") + :documentation ; FDL + "The inverse cumulative distribution function Q(x) for the Gaussian + distribution with standard deviation sigma.") (defun-gsl ugaussian-P (x) "gsl_cdf_ugaussian_P" ((x :double)) :c-return :double - :documentation - "The cumulative distribution function @math{P(x)} for the Gaussian + :documentation ; FDL + "The cumulative distribution function P(x) for the Gaussian distribution with unit standard deviation.") (defun-gsl ugaussian-Q (x) "gsl_cdf_ugaussian_Q" ((x :double)) :c-return :double - :documentation - "The cumulative distribution function @math{Q(x)} for the Gaussian + :documentation ; FDL + "The cumulative distribution function Q(x) for the Gaussian distribution with unit standard deviation.") (defun-gsl ugaussian-Pinv (P) "gsl_cdf_ugaussian_Pinv" ((P :double)) :c-return :double - :documentation - "The inverse cumulative distribution function @math{P(x)} for the Gaussian + :documentation ; FDL + "The inverse cumulative distribution function P(x) for the Gaussian distribution with unit standard deviation.") (defun-gsl ugaussian-Qinv (Q) "gsl_cdf_ugaussian_Qinv" ((Q :double)) :c-return :double - :documentation - "The inverse cumulative distribution function @math{Q(x)} for the Gaussian + :documentation ; FDL + "The inverse cumulative distribution function Q(x) for the Gaussian distribution with unit standard deviation.") ;;; Examples and unit test @@ -137,11 +133,10 @@ "-0.239671528273d+02" "-0.679280164729d+01" "-0.390913184336d+00" "0.893555545521d+01" "-0.176477945898d+00") (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 - collect - (gaussian *rng-mt19937* 10.0d0))))) + collect + (gaussian rng 10.0d0))))) (lisp-unit:assert-first-fp-equal "0.398942280401d-01" (gaussian-pdf 0.0d0 10.0d0)) @@ -151,11 +146,10 @@ "-0.143332423262d+02" "0.945508871923d+00" "-0.111203439859d+02" "-0.307919574280d+01" "-0.209068077307d+01") (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 - collect - (gaussian-ziggurat *rng-mt19937* 10.0d0))))) + collect + (gaussian-ziggurat rng 10.0d0))))) ;; Given in examples in GSL documentation (lisp-unit:assert-first-fp-equal "0.977249868052d+00" (ugaussian-p 2.0d0)) (lisp-unit:assert-first-fp-equal "0.227501319482d-01" (ugaussian-q 2.0d0)) diff --git a/random/generators.lisp b/random/generators.lisp index 7b46f30b60e02631c5731f49fd5d2e2241939e24..52f412dbfd2507e57017783fd6a66d932fbbe535 100644 --- a/random/generators.lisp +++ b/random/generators.lisp @@ -1,6 +1,6 @@ ;; Generators of random numbers. ;; Liam Healy, Sat Jul 15 2006 - 14:43 -;; Time-stamp: <2008-02-02 23:01:02EST generators.lisp> +;; Time-stamp: <2008-02-03 11:36:39EST generators.lisp> ;; $Id: $ (in-package :gsl) @@ -244,10 +244,6 @@ ;;;; Examples and unit test ;;;;**************************************************************************** -;;; (defparameter *rng-mt19937* (make-random-number-generator *mt19937*)) -;;; (defparameter *rng-cmrg* (make-random-number-generator *cmrg*)) -;;; (defparameter *rng-default* (make-random-number-generator *default-type*)) - (lisp-unit:define-test random-number-generators (lisp-unit:assert-equal '(999 162 282 947 231 484 957 744 540 739 759) diff --git a/random/geometric.lisp b/random/geometric.lisp index b94dbac3eaabb2b47221750646a24015416f1140..fde5be5970cafabdfd5dd892346e1fb202a90aeb 100644 --- a/random/geometric.lisp +++ b/random/geometric.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: geometric.lisp -; description: Geometric distribution -; date: Sat Nov 25 2006 - 16:00 -; author: Liam M. Healy -; modified: Sat Nov 25 2006 - 21:56 -;******************************************************** -;;; $Id: $ +;; Geometric distribution +;; Liam Healy, Sat Nov 25 2006 - 16:00 +;; Time-stamp: <2008-02-03 11:04:00EST geometric.lisp> +;; $Id: $ (in-package :gsl) @@ -13,44 +9,45 @@ "gsl_ran_geometric" (((generator generator) :pointer) (p :double)) :c-return :uint - :documentation + :documentation ; FDL "A random integer from the geometric distribution, - the number of independent trials with probability @var{p} until the + the number of independent trials with probability p until the first success. The probability distribution for geometric variates - is p(k) = p (1-p)^{k-1} for @math{k >= 1}. - Note that the distribution begins with @math{k=1} with this - definition. There is another convention in which the exponent @math{k-1} - is replaced by @math{k}.") + is p(k) = p (1-p)^{k-1} for k >= 1. + Note that the distribution begins with k=1 with this + definition. There is another convention in which the exponent k-1 + is replaced by k.") (defun-gsl geometric-pdf (k p) "gsl_ran_geometric_pdf" ((k :uint) (p :double)) :c-return :double - :documentation - "The probability @math{p(k)} of obtaining @var{k} - from a geometric distribution with probability parameter @var{p}, using + :documentation ; FDL + "The probability p(k) of obtaining k + from a geometric distribution with probability parameter p, using the formula given in #'geometric.") (defun-gsl geometric-P (k p) "gsl_cdf_geometric_P" ((k :uint) (p :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{P(k)} for the Geometric distribution with parameter @var{p}.") + :documentation ; FDL + "The cumulative distribution functions + P(k) for the geometric distribution with parameter p.") (defun-gsl geometric-Q (k p) "gsl_cdf_geometric_Q" ((k :uint) (p :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{Q(k)} for the Geometric distribution with parameters @var{p}.") + :documentation ; FDL + "The cumulative distribution functions + Q(k) for the geometric distribution with parameters p.") ;;; Examples and unit test (lisp-unit:define-test geometric (lisp-unit:assert-equal '(1 4 3 1 3 2 1 1 2 1 1) - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 - collect - (geometric *rng-mt19937* 0.4d0)))) + collect + (geometric rng 0.4d0)))) (lisp-unit:assert-first-fp-equal "0.240000000000d+00" (geometric-pdf 2 0.4d0)) diff --git a/random/gumbel1.lisp b/random/gumbel1.lisp index d259d6f090d49845fc3425bf00a7bcc13995eeae..61f207ac27e320d9f12e160abe165500ed3d3357 100644 --- a/random/gumbel1.lisp +++ b/random/gumbel1.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: gumbel1.lisp -; description: Beta distribution -; date: Sun Oct 29 2006 -; author: Liam M. Healy -; modified: Sat Nov 11 2006 - 20:58 -;******************************************************** -;;; $Id: $ +;; The Gumbel type 1 random number distribution +;; Liam Healy, Sun Oct 29 2006 +;; Time-stamp: <2008-02-03 09:46:20EST gumbel1.lisp> +;; $Id: $ (in-package :gsl) @@ -13,47 +9,51 @@ "gsl_ran_gumbel1" (((generator generator) :pointer) (a :double) (b :double)) :c-return :double - :documentation + :documentation ; FDL "A random variate from the Type-1 Gumbel distribution, p(x) dx = a b \exp(-(b \exp(-ax) + ax)) dx - for @math{-\infty < x < \infty}.") + for -\infty < x < \infty.") (defun-gsl gumbel1-pdf (x a b) "gsl_ran_gumbel1_pdf" ((x :double) (a :double) (b :double)) :c-return :double - :documentation - "The probability density @math{p(x)} at @var{x} - for a Type-1 Gumbel distribution with parameters @var{a} and @var{b}, + :documentation ; FDL + "The probability density p(x) at x + for a Type-1 Gumbel distribution with parameters a and b, using the formula given for #'gumbel1.") (defun-gsl gumbel1-P (x a b) "gsl_cdf_gumbel1_P" ((x :double) (a :double) (b :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{P(x)} for the Type-1 Gumbel distribution with - parameters @var{a} and @var{b}.") + :documentation ; FDL + "The cumulative distribution functions + P(x) for the Type-1 Gumbel distribution with + parameters a and b.") (defun-gsl gumbel1-Q (x a b) "gsl_cdf_gumbel1_Q" ((x :double) (a :double) (b :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{Q(x)} for the Type-1 Gumbel distribution with - parameters @var{a} and @var{b}.") + :documentation ; FDL + "The cumulative distribution functions + Q(x) for the Type-1 Gumbel distribution with + parameters a and b.") (defun-gsl gumbel1-Pinv (P a b) "gsl_cdf_gumbel1_Pinv" ((P :double) (a :double) (b :double)) :c-return :double - :documentation "The inverse cumulative distribution functions - @math{P(x)} for the Type-1 Gumbel distribution with - parameters @var{a} and @var{b}.") + :documentation ; FDL + "The inverse cumulative distribution functions + P(x) for the Type-1 Gumbel distribution with + parameters a and b.") (defun-gsl gumbel1-Qinv (Q a b) "gsl_cdf_gumbel1_Qinv" ((Q :double) (a :double) (b :double)) :c-return :double - :documentation "The inverse cumulative distribution functions - @math{Q(x)} for the Type-1 Gumbel distribution with - parameters @var{a} and @var{b}.") + :documentation ; FDL + "The inverse cumulative distribution functions + Q(x) for the Type-1 Gumbel distribution with + parameters a and b.") ;;; Examples and unit test (lisp-unit:define-test gumbel1 diff --git a/random/gumbel2.lisp b/random/gumbel2.lisp index 7a4be3cdf3eceffae28d54987ab6f8d731c08eaf..23109738f1484631413666a510b2714a0187d80d 100644 --- a/random/gumbel2.lisp +++ b/random/gumbel2.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: gumbel2.lisp -; description: Beta distribution -; date: Sun Oct 29 2006 -; author: Liam M. Healy -; modified: Sat Nov 11 2006 - 20:58 -;******************************************************** -;;; $Id: $ +;; The Gumbel type 2 random number distribution +;; Liam Healy, Sun Oct 29 2006 +;; Time-stamp: <2008-02-03 10:10:26EST gumbel2.lisp> +;; $Id: $ (in-package :gsl) @@ -13,47 +9,49 @@ "gsl_ran_gumbel2" (((generator generator) :pointer) (a :double) (b :double)) :c-return :double - :documentation + :documentation ; FDL "A random variate from the Type-2 Gumbel - distribution, - p(x) dx = a b x^@{-a-1@} \exp(-b x^@{-a@}) dx - for @math{0 < x < \infty}.") + distribution, p(x) dx = a b x^{-a-1} \exp(-b x^{-a}) dx + for 0 < x < \infty.") (defun-gsl gumbel2-pdf (x a b) "gsl_ran_gumbel2_pdf" ((x :double) (a :double) (b :double)) :c-return :double - :documentation - "The probability density @math{p(x)} at @var{x} - for a Type-2 Gumbel distribution with parameters @var{a} and @var{b}, + :documentation ; FDL + "The probability density p(x) at x + for a Type-2 Gumbel distribution with parameters a and b, using the formula given in #'gumbel2.") (defun-gsl gumbel2-P (x a b) "gsl_cdf_gumbel2_P" ((x :double) (a :double) (b :double)) - :c-return :double + :c-return :double ; FDL :documentation "The cumulative distribution functions - @math{P(x)} for the Type-2 Gumbel distribution with - parameters @var{a} and @var{b}.") + P(x) for the Type-2 Gumbel distribution with + parameters a and b.") (defun-gsl gumbel2-Q (x a b) "gsl_cdf_gumbel2_Q" ((x :double) (a :double) (b :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{Q(x)} for the Type-2 Gumbel distribution with - parameters @var{a} and @var{b}.") + :documentation ; FDL + "The cumulative distribution functions + Q(x) for the Type-2 Gumbel distribution with + parameters a and b.") (defun-gsl gumbel2-Pinv (P a b) "gsl_cdf_gumbel2_Pinv" ((P :double) (a :double) (b :double)) :c-return :double - :documentation "The inverse cumulative distribution functions - @math{P(x)} for the Type-2 Gumbel distribution with - parameters @var{a} and @var{b}.") + :documentation ; FDL + "The inverse cumulative distribution functions + P(x) for the Type-2 Gumbel distribution with + parameters a and b.") (defun-gsl gumbel2-Qinv (Q a b) "gsl_cdf_gumbel2_Qinv" ((Q :double) (a :double) (b :double)) :c-return :double - :documentation "The inverse cumulative distribution functions - @math{Q(x)} for the Type-2 Gumbel distribution with - parameters @var{a} and @var{b}.") + :documentation ; FDL + "The inverse cumulative distribution functions + Q(x) for the Type-2 Gumbel distribution with + parameters a and b.") ;;; Examples and unit test (lisp-unit:define-test gumbel2 @@ -63,11 +61,10 @@ "0.460260800603d+02" "0.677268352507d+01" "0.324619836866d+01" "0.664079780729d+01" "0.728568789702d+01") (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 - collect - (gumbel2 *rng-mt19937* 1.0d0 2.0d0))))) + collect + (gumbel2 rng 1.0d0 2.0d0))))) (lisp-unit:assert-first-fp-equal "0.536256036829d-01" (gumbel2-pdf 5.0d0 1.0d0 2.0d0)) diff --git a/random/hypergeometric.lisp b/random/hypergeometric.lisp index fc1aa1b05fe69e8f585486f9af4bdf6fb995c3c2..242c48b707d1acdd9e204458dc3a66435ebac8ad 100644 --- a/random/hypergeometric.lisp +++ b/random/hypergeometric.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: hypergeometric.lisp -; description: Hypergeometric distribution -; date: Sat Nov 25 2006 - 16:00 -; author: Liam M. Healy -; modified: Sat Nov 25 2006 - 22:40 -;******************************************************** -;;; $Id: $ +;; Hypergeometric distribution +;; Liam Healy, Sat Nov 25 2006 - 16:00 +;; Time-stamp: <2008-02-03 10:12:46EST hypergeometric.lisp> +;; $Id: $ (in-package :gsl) @@ -13,53 +9,50 @@ "gsl_ran_hypergeometric" (((generator generator) :pointer) (n1 :uint) (n2 :uint)(tt :uint)) :c-return :uint - :documentation + :documentation ; FDL "A random integer from the hypergeometric distribution. The probability distribution for hypergeometric random variates is p(k) = C(n_1, k) C(n_2, t - k) / C(n_1 + n_2, t) - where @math{C(a,b) = a!/(b!(a-b)!)} and - @math{t <= n_1 + n_2}. The domain of @math{k} is - @math{max(0,t-n_2), ..., min(t,n_1)}. - If a population contains @math{n_1} elements of ``type 1'' and - @math{n_2} elements of ``type 2'' then the hypergeometric - distribution gives the probability of obtaining @math{k} elements of - ``type 1'' in @math{t} samples from the population without + where C(a,b) = a!/(b!(a-b)!) and + t <= n_1 + n_2. The domain of k is + max(0,t-n_2), ..., min(t,n_1). + If a population contains n_1 elements of ``type 1'' and + n_2 elements of ``type 2'' then the hypergeometric + distribution gives the probability of obtaining k elements of + ``type 1'' in t samples from the population without replacement.") (defun-gsl hypergeometric-pdf (k n1 n2 tt) "gsl_ran_hypergeometric_pdf" ((k :uint) (n1 :uint) (n2 :uint)(tt :uint)) :c-return :double - :documentation - "The probability @math{p(k)} of obtaining @var{k} - from a hypergeometric distribution with parameters @var{n1}, @var{n2}, - @var{tt}, using the formula given in #'hypergeometric.") + :documentation ; FDL + "The probability p(k) of obtaining k + from a hypergeometric distribution with parameters n1, n2, + tt, using the formula given in #'hypergeometric.") (defun-gsl hypergeometric-P (k n1 n2 tt) "gsl_cdf_hypergeometric_P" ((k :uint) (n1 :uint) (n2 :uint)(tt :uint)) :c-return :double - :documentation - "The cumulative distribution functions @math{P(k)} for the - hypergeometric distribution with parameters @var{n1}, @var{n2} - and @var{tt}.") + :documentation ; FDL + "The cumulative distribution functions P(k) for the + hypergeometric distribution with parameters n1, n2 and tt.") (defun-gsl hypergeometric-Q (k n1 n2 tt) "gsl_cdf_hypergeometric_Q" ((k :uint) (n1 :uint) (n2 :uint)(tt :uint)) :c-return :double - :documentation - "The cumulative distribution functions @math{Q(k)} for the - hypergeometric distribution with parameters @var{n1}, @var{n2} - and @var{tt}.") + :documentation ; FDL + "The cumulative distribution functions Q(k) for the + hypergeometric distribution with parameters n1, n2, and tt.") ;;; Examples and unit test (lisp-unit:define-test hypergeometric-randist (lisp-unit:assert-equal '(2 1 0 0 1 1 3 1 0 1 3) - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 - collect - (hypergeometric *rng-mt19937* 3 6 3)))) + collect + (hypergeometric rng 3 6 3)))) (lisp-unit:assert-first-fp-equal "0.357142857143d+00" (hypergeometric-pdf 0 2 6 3)) diff --git a/random/landau.lisp b/random/landau.lisp index 9d63ffa6f528987f9b66be69b4ffb5f8e4cfb03a..d15104b15e6f16c472dbaffc62ada8aed535d0b2 100644 --- a/random/landau.lisp +++ b/random/landau.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: landau.lisp -; description: Landau distribution -; date: Sat Sep 30 2006 -; author: Liam M. Healy -; modified: Sat Sep 30 2006 - 21:52 -;******************************************************** -;;; $Id: $ +;; Landau distribution +;; Liam Healy, Sat Sep 30 2006 +;; Time-stamp: <2008-02-03 11:30:46EST landau.lisp> +;; $Id: $ (in-package :gsl) @@ -13,7 +9,7 @@ "gsl_ran_landau" (((generator generator) :pointer)) :c-return :double - :documentation + :documentation ; FDL "A random variate from the Landau distribution. The probability distribution for Landau random variates is defined analytically by the complex integral, @@ -25,8 +21,8 @@ (defun-gsl landau-pdf (x) "gsl_ran_landau_pdf" ((x :double)) :c-return :double - :documentation - "The probability density @math{p(x)} at @var{x} + :documentation ; FDL + "The probability density p(x) at x for the Landau distribution using an approximation to the formula given in #'landau.") @@ -38,11 +34,9 @@ "0.261461684796d+02" "0.433721764097d+01" "0.167995462811d+01" "0.424757192183d+01" "0.468150620898d+01") (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 - collect - (landau *rng-mt19937*))))) + collect (landau rng))))) (lisp-unit:assert-first-fp-equal "0.173319689959d+00" (landau-pdf 0.25d0))) diff --git a/random/laplace.lisp b/random/laplace.lisp index 778b9893bb0612d1199e741f981c0b7bf37c47a2..343ddf320ce91a6f5d8956d4998db00ee609c43e 100644 --- a/random/laplace.lisp +++ b/random/laplace.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: exponential.lisp -; description: Exponential distribution -; date: Sun Sep 17 2006 -; author: Liam M. Healy -; modified: Sat Sep 30 2006 - 19:17 -;******************************************************** -;;; $Id: $ +;; Exponential distribution +;; Liam Healy, Sun Sep 17 2006 +;; Time-stamp: <2008-02-03 11:06:31EST laplace.lisp> +;; $Id: $ (in-package :gsl) @@ -13,43 +9,47 @@ "gsl_ran_laplace" (((generator generator) :pointer) (a :double)) :c-return :double - :documentation - "A random variate from the Laplace distribution with width @var{a}. + :documentation ; FDL + "A random variate from the Laplace distribution with width a. The distribution is p(x) dx = {1 \over 2 a} \exp(-|x/a|) dx - for @math{-\infty < x < \infty}.") + for -\infty < x < \infty.") (defun-gsl laplace-pdf (x a) "gsl_ran_laplace_pdf" ((x :double) (a :double)) :c-return :double - :documentation - "The probability density @math{p(x)} at @var{x} - for a Laplace distribution with width @var{a}, using the formula + :documentation ; FDL + "The probability density p(x) at x + for a Laplace distribution with width a, using the formula given for #'laplace.") (defun-gsl laplace-P (x a) "gsl_cdf_laplace_P" ((x :double) (a :double)) :c-return :double - :documentation "The cumulative distribution function - @math{P(x)} for the laplace distribution with width @var{a}.") + :documentation ; FDL + "The cumulative distribution function + P(x) for the laplace distribution with width a.") (defun-gsl laplace-Q (x a) "gsl_cdf_laplace_Q" ((x :double) (a :double)) :c-return :double - :documentation "The cumulative distribution function - @math{Q(x)} for the laplace distribution with width @var{a}.") + :documentation ; FDL + "The cumulative distribution function + Q(x) for the laplace distribution with width a.") (defun-gsl laplace-Pinv (P a) "gsl_cdf_laplace_Pinv" ((P :double) (a :double)) :c-return :double - :documentation "The inverse cumulative distribution function - @math{P(x)} for the laplace distribution with width @var{a}.") + :documentation ; FDL + "The inverse cumulative distribution function + P(x) for the laplace distribution with width a.") (defun-gsl laplace-Qinv (Q a) "gsl_cdf_laplace_Qinv" ((Q :double) (a :double)) :c-return :double - :documentation "The inverse cumulative distribution function - @math{Q(x)} for the laplace distribution with width @var{a}.") + :documentation ; FDL + "The inverse cumulative distribution function + Q(x) for the laplace distribution with width a.") ;;; Examples and unit test (lisp-unit:define-test laplace @@ -59,11 +59,10 @@ "0.888815832003d+00" "0.716189249197d+01" "0.252463778091d+02" "0.734165104806d+01" "0.654142651602d+01") (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 collect - (laplace *rng-mt19937* 10.0d0))))) + (laplace rng 10.0d0))))) (lisp-unit:assert-first-fp-equal "0.500000000000d-01" (laplace-pdf 0.0d0 10.0d0)) diff --git a/random/levy.lisp b/random/levy.lisp index 8d5eb22f9d2723c74fae082a9fac3d29eae2d245..969887b19e578c2a6a63ad85ff77f91511a71ac8 100644 --- a/random/levy.lisp +++ b/random/levy.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: levy.lisp -; description: Levy distribution -; date: Sat Sep 30 2006 -; author: Liam M. Healy -; modified: Sat Sep 30 2006 - 22:45 -;******************************************************** -;;; $Id: $ +;; Levy distribution +;; Liam Healy, Sat Sep 30 2006 +;; Time-stamp: <2008-02-03 10:44:49EST levy.lisp> +;; $Id: $ (in-package :gsl) @@ -13,39 +9,38 @@ "gsl_ran_levy" (((generator generator) :pointer) (c :double) (alpha :double)) :c-return :double - :documentation + :documentation ; FDL "A random variate from the Levy symmetric stable - distribution with scale @var{c} and exponent @var{alpha}. The symmetric - stable probability distribution is defined by a fourier transform + distribution with scale c and exponent alpha. The symmetric + stable probability distribution is defined by a fourier transform, p(x) = {1 \over 2 \pi} \int_{-\infty}^{+\infty} dt \exp(-it x - |c t|^\alpha) - There is no explicit solution for the form of @math{p(x)} and the - library does not define a corresponding @code{pdf} function. For - @math{\alpha = 1} the distribution reduces to the Cauchy distribution. For - @math{\alpha = 2} it is a Gaussian distribution with @c{$\sigma = \sqrt{2} c$} - @math{\sigma = \sqrt@{2@} c}. For @math{\alpha < 1} the tails of the - distribution become extremely wide. - The algorithm only works for @c{$0 < \alpha \le 2$} @math{0 < alpha <= 2}.") + There is no explicit solution for the form of p(x) and the + library does not define a corresponding pdf function. For + \alpha = 1 the distribution reduces to the Cauchy distribution. For + \alpha = 2 it is a Gaussian distribution with \sigma = \sqrt{2} c + For \alpha < 1 the tails of the distribution become extremely wide. + The algorithm only works for 0 < alpha <= 2.") (defun-gsl levy-skew (generator c alpha beta) "gsl_ran_levy_skew" (((generator generator) :pointer) (c :double) (alpha :double) (beta :double)) :c-return :double - :documentation + :documentation ; FDL "A random variate from the Levy skew stable - distribution with scale @var{c}, exponent @var{alpha} and skewness - parameter @var{beta}. The skewness parameter must lie in the range - @math{[-1,1]}. The Levy skew stable probability distribution is defined + distribution with scale c exponent alpha and skewness + parameter beta. The skewness parameter must lie in the range + [-1,1]. The Levy skew stable probability distribution is defined by a fourier transform, - p(x) = {1 \over 2 \pi} \int_{-\infty}^{+\infty} dt \exp(-it x - |c t|^\alpha (1-i \beta \sign(t) \tan(\pi\alpha/2))) - When @math{\alpha = 1} the term @math{\tan(\pi \alpha/2)} is replaced by - @math{-(2/\pi)\log|t|}. There is no explicit solution for the form of - @math{p(x)} and the library does not define a corresponding @code{pdf} - function. For @math{\alpha = 2} the distribution reduces to a Gaussian - distribution with @c{$\sigma = \sqrt{2} c$} - @math{\sigma = \sqrt@{2@} c} and the skewness parameter has no effect. - For @math{\alpha < 1} the tails of the distribution become extremely - wide. The symmetric distribution corresponds to @math{\beta = 0}. - The algorithm only works for @c{$0 < \alpha \le 2$} @math{0 < alpha <= 2}.") + p(x) = {1 \over 2 \pi} \int_{-\infty}^{+\infty} dt + \exp(-it x - |c t|^\alpha (1-i \beta \sign(t) \tan(\pi\alpha/2))) + When \alpha = 1 the term \tan(\pi \alpha/2) is replaced by + -(2/\pi)\log|t|. There is no explicit solution for the form of + p(x)} and the library does not define a corresponding pdf + function. For \alpha = 2 the distribution reduces to a Gaussian + distribution with \sigma = \sqrt{2} c and the skewness parameter + has no effect. For \alpha < 1 the tails of the distribution + become extremely wide. The symmetric distribution corresponds to \beta = 0. + The algorithm only works for 0 < \alpha \le 2.") ;;; Examples and unit test @@ -56,19 +51,17 @@ "-0.510713467479d+00" "0.164820785369d+00" "-0.148578990410d+00" "-0.190748857444d+01" "-0.208619521400d+01") (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 collect - (levy *rng-mt19937* 1.0d0 2.0d0))))) + (levy rng 1.0d0 2.0d0))))) (lisp-unit:assert-equal '("0.269410983324d+01" "-0.293954386447d+00" "-0.127034013523d+01" "0.107715386401d+01" "0.137712184069d+00" "0.941972843811d+00" "-0.510713467479d+00" "0.164820785369d+00" "-0.148578990410d+00" "-0.190748857444d+01" "-0.208619521400d+01") (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 collect - (levy-skew *rng-mt19937* 1.0d0 2.0d0 1.0d0)))))) + (levy-skew rng 1.0d0 2.0d0 1.0d0)))))) diff --git a/random/logarithmic.lisp b/random/logarithmic.lisp index 8a7228024b2d95fde636ff6c31b553fab254456a..e034879a5a5f749bc8b1c375dd552d8fef96dcd2 100644 --- a/random/logarithmic.lisp +++ b/random/logarithmic.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: logarithmic.lisp -; description: Logarithmic distribution -; date: Sat Nov 25 2006 - 16:00 -; author: Liam M. Healy -; modified: Sat Nov 25 2006 - 22:57 -;******************************************************** -;;; $Id: $ +;; Logarithmic distribution +;; Liam Healy, Sat Nov 25 2006 - 16:00 +;; Time-stamp: <2008-02-03 11:10:27EST logarithmic.lisp> +;; $Id: $ (in-package :gsl) @@ -13,29 +9,28 @@ "gsl_ran_logarithmic" (((generator generator) :pointer) (p :double)) :c-return :uint - :documentation + :documentation ; FDL "A random integer from the logarithmic distribution. The probability distribution for logarithmic random variates is p(k) = {-1 \over \log(1-p)} {\left( p^k \over k \right)} - for @math{k >= 1}.") + for k >= 1.") (defun-gsl logarithmic-pdf (k p) "gsl_ran_logarithmic_pdf" ((k :uint) (p :double)) :c-return :double :documentation - "The probability @math{p(k)} of obtaining @var{k} - from a logarithmic distribution with probability parameter @var{p}, + "The probability p(k) of obtaining k + from a logarithmic distribution with probability parameter p, using the formula given in #'logarithmic.") ;;; Examples and unit test (lisp-unit:define-test logarithmic (lisp-unit:assert-equal '(1 3 1 4 1 1 2 1 1 5 2) - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 - collect - (logarithmic *rng-mt19937* 0.9d0)))) + collect + (logarithmic rng 0.9d0)))) (lisp-unit:assert-first-fp-equal "0.156609215118d+00" (logarithmic-pdf 2 0.4d0))) diff --git a/random/logistic.lisp b/random/logistic.lisp index 53b31266cb6c7d352a74ab8ee3e18f6f5a49a793..a61b4ddb4bb5cd12745293de4adf94b5655fa239 100644 --- a/random/logistic.lisp +++ b/random/logistic.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: logistic.lisp -; description: logistic distribution -; date: Sat Oct 7 2006 - 16:13 -; author: Liam M. Healy -; modified: Sun Oct 8 2006 - 17:32 -;******************************************************** -;;; $Id: $ +;; Logistic distribution +;; Liam Healy, Sat Oct 7 2006 - 16:13 +;; Time-stamp: <2008-02-03 10:35:08EST logistic.lisp> +;; $Id: $ (in-package :gsl) @@ -13,42 +9,46 @@ "gsl_ran_logistic" (((generator generator) :pointer) (a :double)) :c-return :double - :documentation + :documentation ; FDL "A random variate from the logistic distribution. The distribution function is p(x) dx = { \exp(-x/a) \over a (1 + \exp(-x/a))^2 } dx - for @math{-\infty < x < +\infty}.") + for -\infty < x < +\infty.") (defun-gsl logistic-pdf (x a) "gsl_ran_logistic_pdf" ((x :double) (a :double)) :c-return :double - :documentation - "The probability density @math{p(x)} at @var{x} - for a logistic distribution with scale parameter @var{a}, using the + :documentation ; FDL + "The probability density p(x) at x + for a logistic distribution with scale parameter a, using the formula given in #'logistic.") (defun-gsl logistic-P (x a) "gsl_cdf_logistic_P" ((x :double) (a :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{P(x)} for the logistic distribution with scale parameter @var{a}.") + :documentation ; FDL + "The cumulative distribution functions + P(x) for the logistic distribution with scale parameter a.") (defun-gsl logistic-Q (x a) "gsl_cdf_logistic_Q" ((x :double) (a :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{Q(x)} for the logistic distribution with scale parameter @var{a}.") + :documentation ; FDL + "The cumulative distribution functions + Q(x) for the logistic distribution with scale parameter a.") (defun-gsl logistic-Pinv (P a) "gsl_cdf_logistic_Pinv" ((P :double) (a :double)) :c-return :double - :documentation "The inverse cumulative distribution functions - @math{P(x)} for the logistic distribution with scale parameter @var{a}.") + :documentation ; FDL + "The inverse cumulative distribution functions + P(x) for the logistic distribution with scale parameter a.") (defun-gsl logistic-Qinv (Q a) "gsl_cdf_logistic_Qinv" ((Q :double) (a :double)) :c-return :double - :documentation "The inverse cumulative distribution functions - @math{Q(x)} for the logistic distribution with scale parameter @var{a}.") + :documentation ; FDL + "The inverse cumulative distribution functions + Q(x) for the logistic distribution with scale parameter a.") ;;; Examples and unit test (lisp-unit:define-test logistic @@ -58,11 +58,10 @@ "0.311425552636d+02" "0.106846737210d+02" "0.160518409540d+01" "0.104572419047d+02" "0.115237141063d+02") (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 collect - (logistic *rng-mt19937* 10.0d0))))) + (logistic rng 10.0d0))))) (lisp-unit:assert-first-fp-equal "0.235003712202d+00" (logistic-pdf 0.5d0 1.0d0)) diff --git a/random/lognormal.lisp b/random/lognormal.lisp index b23c979641ad327c5478c53d599dea7c81771067..a0dd11702110e58dabc5363e03ec688cb93b5b7b 100644 --- a/random/lognormal.lisp +++ b/random/lognormal.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: lognormal.lisp -; description: Lognormal distribution -; date: Sat Sep 30 2006 -; author: Liam M. Healy -; modified: Sat Oct 7 2006 - 15:18 -;******************************************************** -;;; $Id: $ +;; Lognormal distribution +;; Liam Healy, Sat Sep 30 2006 +;; Time-stamp: <2008-02-03 11:08:34EST lognormal.lisp> +;; $Id: $ (in-package :gsl) @@ -13,47 +9,50 @@ "gsl_ran_lognormal" (((generator generator) :pointer) (zeta :double) (sigma :double)) :c-return :double - :documentation + :documentation ; FDL "A random variate from the lognormal distribution. The distribution function is p(x) dx = {1 \over x \sqrt{2 \pi \sigma^2}} \exp(-(\ln(x) - \zeta)^2/2 \sigma^2) dx - for @math{x > 0}.") + for x > 0.") (defun-gsl lognormal-pdf (x zeta sigma) "gsl_ran_lognormal_pdf" ((x :double) (zeta :double) (sigma :double)) :c-return :double - :documentation - "The probability density @math{p(x)} at @var{x} - for a lognormal distribution with parameters @var{zeta} and @var{sigma}, + :documentation ; FDL + "The probability density p(x) at X + for a lognormal distribution with parameters zeta and sigma, using the formula given in #'lognormal.") (defun-gsl lognormal-P (x zeta sigma) "gsl_cdf_lognormal_P" ((x :double) (zeta :double) (sigma :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{P(x)} for the lognormal distribution with parameters - @var{zeta} and @var{sigma}.") + :documentation ; FDL + "The cumulative distribution functions + P(x) for the lognormal distribution with parameters zeta and sigma.") (defun-gsl lognormal-Q (x zeta sigma) "gsl_cdf_lognormal_Q" ((x :double) (zeta :double) (sigma :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{Q(x)} for the lognormal distribution with parameters - @var{zeta} and @var{sigma}.") + :documentation ; FDL + "The cumulative distribution functions + Q(x) for the lognormal distribution with parameters + zeta and sigma.") (defun-gsl lognormal-Pinv (P zeta sigma) "gsl_cdf_lognormal_Pinv" ((P :double) (zeta :double) (sigma :double)) :c-return :double - :documentation "The inverse cumulative distribution functions - @math{P(x)} for the lognormal distribution with parameters - @var{zeta} and @var{sigma}.") + :documentation ; FDL + "The inverse cumulative distribution functions + P(x) for the lognormal distribution with parameters + zeta and sigma.") (defun-gsl lognormal-Qinv (Q zeta sigma) "gsl_cdf_lognormal_Qinv" ((Q :double) (zeta :double) (sigma :double)) :c-return :double - :documentation "The inverse cumulative distribution functions - @math{Q(x)} for the lognormal distribution with parameters - @var{zeta} and @var{sigma}.") + :documentation ; FDL + "The inverse cumulative distribution functions + Q(x) for the lognormal distribution with parameters + zeta and sigma.") ;;; Examples and unit test (lisp-unit:define-test lognormal @@ -63,11 +62,10 @@ "0.914620656772d+00" "0.727307901066d+00" "0.218018485218d+01" "0.389088566169d+01" "0.182184697889d+03") (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 collect - (lognormal *rng-mt19937* 1.0d0 2.0d0))))) + (lognormal rng 1.0d0 2.0d0))))) (lisp-unit:assert-first-fp-equal "0.152898339657d+00" (lognormal-pdf 1.2d0 1.0d0 2.0d0)) diff --git a/random/multinomial.lisp b/random/multinomial.lisp index 37ce552e1be7dfa57cb0bf552ad1c3d8829cb422..68c5c34451506cdae28aa5eff8a9362568ae4e6d 100644 --- a/random/multinomial.lisp +++ b/random/multinomial.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: multinomial.lisp -; description: Multinomial distribution -; date: Sat Nov 25 2006 - 16:00 -; author: Liam M. Healy -; modified: Sat Nov 25 2006 - 21:16 -;******************************************************** -;;; $Id: $ +;; Multinomial distribution +;; Liam Healy, Sat Nov 25 2006 - 16:00 +;; Time-stamp: <2008-02-03 10:28:38EST multinomial.lisp> +;; $Id: $ (in-package :gsl) @@ -18,63 +14,57 @@ ;; technically, n should be a uint array, but integers work ((gsl-array n) :pointer)) :c-return :void - :documentation - "Returns an array n of @var{K} random variates from a + :documentation ; FDL + "Returns an array n of K random variates from a multinomial distribution. The sum of the array n is specified by sum=N. The distribution function is P(n_1, n_2, ..., n_K) = (N!/(n_1! n_2! ... n_K!)) p_1^n_1 p_2^n_2 ... p_K^n_K - where @math{(n_1, n_2, ..., n_K)} - are nonnegative integers with - @math{sum_@{k=1@}^K n_k = N}, - and @math{(p_1, p_2, ..., p_K)} - is a probability distribution with @math{\sum p_i = 1}. - If the array @var{p}[@var{K}] is not normalized then its entries will be + where (n_1, n_2, ..., n_K) are nonnegative integers with + sum_@{k=1@}^K n_k = N, and (p_1, p_2, ..., p_K) + is a probability distribution with \sum p_i = 1. + If the array p[K] is not normalized then its entries will be treated as weights and normalized appropriately. Random variates are generated using the conditional binomial method (see - C.S. David, @cite{The computer generation of multinomial random - variates}, Comp. Stat. Data Anal. 16 (1993) 205--217 for details).") + C.S. David, \"The computer generation of multinomial random + variates,\" Comp. Stat. Data Anal. 16 (1993) 205--217 for details).") (defun-gsl multinomial-pdf (p n) "gsl_ran_multinomial_pdf" (((dim0 p) :uint) ((gsl-array p) :pointer) ((gsl-array n) :pointer)) :c-return :double - :documentation - "Compute the probability @math{P(n_1, n_2, ..., n_K)} - of sampling @var{n}[@var{K}] from a multinomial distribution - with parameters @var{p}[@var{K}], using the formula given for #'multinomial.") + :documentation ; FDL + "Compute the probability P(n_1, n_2, ..., n_K) + of sampling n[K] from a multinomial distribution + with parameters p[K], using the formula given for #'multinomial.") (defun-gsl multinomial-log-pdf (p n) "gsl_ran_multinomial_lnpdf" (((dim0 p) :uint) ((gsl-array p) :pointer) ((gsl-array n) :pointer)) :c-return :double - :documentation - "Compute the natural logarithm of the probability @math{P(n_1, n_2, ..., n_K)} - of sampling @var{n}[@var{K}] from a multinomial distribution - with parameters @var{p}[@var{K}], using the formula given for #'multinomial.") + :documentation ; FDL + "Compute the natural logarithm of the probability P(n_1, n_2, ..., n_K) + of sampling n[K] from a multinomial distribution + with parameters p[K], using the formula given for #'multinomial.") ;;; Examples and unit test (lisp-unit:define-test multinomial (lisp-unit:assert-equalp #(5 0 1 2) - (progn - (rng-set *rng-mt19937* 0) - (with-data (p vector-double 4) - (setf (data p) #(0.1d0 0.2d0 0.3d0 0.4d0)) - (with-data (n vector-fixnum 4) - (multinomial *rng-mt19937* 8 p n) - (data n))))) + (letm ((rng (random-number-generator *mt19937* 0)) + (p (vector-double #(0.1d0 0.2d0 0.3d0 0.4d0))) + (n (vector-fixnum 4))) + (multinomial rng 8 p n) + (data n))) (lisp-unit:assert-first-fp-equal "0.806400000000d-04" - (with-data (p vector-double 4) - (setf (data p) #(0.1d0 0.2d0 0.3d0 0.4d0)) - (with-data (n vector-fixnum 4) - (setf (data n) #(5 0 1 2)) - (multinomial-pdf p N)))) + (letm ((p (vector-double #(0.1d0 0.2d0 0.3d0 0.4d0))) + (n (vector-fixnum 4))) + (setf (data n) #(5 0 1 2)) + (multinomial-pdf p N))) (lisp-unit:assert-first-fp-equal "-0.942551575364d+01" - (with-data (p vector-double 4) - (setf (data p) #(0.1d0 0.2d0 0.3d0 0.4d0)) - (with-data (n vector-fixnum 4) - (setf (data n) #(5 0 1 2)) - (multinomial-log-pdf p n))))) + (letm ((p (vector-double #(0.1d0 0.2d0 0.3d0 0.4d0))) + (n (vector-fixnum 4))) + (setf (data n) #(5 0 1 2)) + (multinomial-log-pdf p n)))) diff --git a/random/negative-binomial.lisp b/random/negative-binomial.lisp index 3fed45e1021569448f01652247a2d510f549949a..1a5b5fb0155f5826a36d27507dd6764412f4a55e 100644 --- a/random/negative-binomial.lisp +++ b/random/negative-binomial.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: negative-binomial.lisp -; description: Negative binomial and Pascal distributions -; date: Sat Nov 25 2006 - 16:00 -; author: Liam M. Healy -; modified: Sat Nov 25 2006 - 21:42 -;******************************************************** -;;; $Id: $ +;; Negative binomial and Pascal distributions +;; Liam Healy, Sat Nov 25 2006 - 16:00 +;; Time-stamp: <2008-02-03 11:25:48EST negative-binomial.lisp> +;; $Id: $ (in-package :gsl) @@ -17,35 +13,37 @@ "gsl_ran_negative_binomial" (((generator generator) :pointer) (p :double) (n :double)) :c-return :uint - :documentation + :documentation ; FDL "A random integer from the negative binomial - distribution, the number of failures occurring before @var{n} successes - in independent trials with probability @var{p} of success. The + distribution, the number of failures occurring before n successes + in independent trials with probability p of success. The probability distribution for negative binomial variates is, p(k) = {\Gamma(n + k) \over \Gamma(k+1) \Gamma(n) } p^n (1-p)^k - Note that @math{n} is not required to be an integer.") + Note that n is not required to be an integer.") (defun-gsl negative-binomial-pdf (k p n) "gsl_ran_negative_binomial_pdf" ((k :uint) (p :double) (n :double)) :c-return :double - :documentation - "The probability @math{p(k)} of obtaining @var{k} - from a negative binomial distribution with parameters @var{p} and - @var{n}, using the formula given in #'negative-binomial.") + :documentation ; FDL + "The probability p(k) of obtaining k + from a negative binomial distribution with parameters p and + n, using the formula given in #'negative-binomial.") (defun-gsl negative-binomial-P (k p n) "gsl_cdf_negative_binomial_P" ((k :uint) (p :double) (n :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{P(k)} for the negative binomial distribution - with parameters @var{p} and @var{n}.") + :documentation ; FDL + "The cumulative distribution functions + P(k) for the negative binomial distribution + with parameters p and n.") (defun-gsl negative-binomial-Q (k p n) "gsl_cdf_negative_binomial_Q" ((k :uint) (p :double) (n :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{Q(k)} for the negative binomial distribution - with parameters @var{p} and @var{n}.") + :documentation ; FDL + "The cumulative distribution functions + Q(k) for the negative binomial distribution + with parameters p and n.") ;;;;**************************************************************************** ;;;; Pascal @@ -55,34 +53,36 @@ "gsl_ran_pascal" (((generator generator) :pointer) (p :double) (n :uint)) :c-return :uint - :documentation + :documentation ; FDL "A random integer from the Pascal distribution. The Pascal distribution is simply a negative binomial distribution with an - integer value of @math{n}. + integer value of n. p(k) = {(n + k - 1)! \over k! (n - 1)! } p^n (1-p)^k - @math{k >= 0}.") + k >= 0.") (defun-gsl pascal-pdf (k p n) "gsl_ran_pascal_pdf" ((k :uint) (p :double) (n :uint)) :c-return :double - :documentation - "The probability @math{p(k)} of obtaining @var{k} - from a Pascal distribution with parameters @var{p} and - @var{n}, using the formula given in #'pascal.") + :documentation ; FDL + "The probability p(k) of obtaining k + from a Pascal distribution with parameters p and + n, using the formula given in #'pascal.") (defun-gsl pascal-P (k p n) "gsl_cdf_pascal_P" ((k :uint) (p :double) (n :uint)) :c-return :double - :documentation "The cumulative distribution functions - @math{P(k)} for the Pascal distribution - with parameters @var{p} and @var{n}.") + :documentation ; FDL + "The cumulative distribution functions + P(k) for the Pascal distribution + with parameters p and n.") (defun-gsl pascal-Q (k p n) "gsl_cdf_pascal_Q" ((k :uint) (p :double) (n :uint)) :c-return :double - :documentation "The cumulative distribution functions - @math{Q(k)} for the Pascal distribution - with parameters @var{p} and @var{n}.") + :documentation ; FDL + "The cumulative distribution functions + Q(k) for the Pascal distribution + with parameters p and n.") ;;;;**************************************************************************** ;;;; Examples and unit test @@ -91,11 +91,10 @@ (lisp-unit:define-test negative-binomial (lisp-unit:assert-equal '(10 7 12 23 20 24 18 12 4 22 15) - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 - collect - (negative-binomial *rng-mt19937* 0.4d0 12.0d0)))) + collect + (negative-binomial rng 0.4d0 12.0d0)))) (lisp-unit:assert-first-fp-equal "0.569847670899d-02" (negative-binomial-pdf 5 0.4d0 12.0d0)) @@ -107,11 +106,10 @@ (negative-binomial-Q 5 0.4d0 12.0d0)) (lisp-unit:assert-equal '(10 7 12 23 20 24 18 12 4 22 15) - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 - collect - (pascal *rng-mt19937* 0.4d0 12)))) + collect + (pascal rng 0.4d0 12)))) (lisp-unit:assert-first-fp-equal "0.569847670899d-02" (pascal-pdf 5 0.4d0 12)) diff --git a/random/pareto.lisp b/random/pareto.lisp index 465b126c6f916a0b8584e359e104e732b9b19f93..f269b1855d4558847a83d4d52233b5ab5f8eb082 100644 --- a/random/pareto.lisp +++ b/random/pareto.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: pareto.lisp -; description: Pareto distribution -; date: Sat Oct 8 2006 - 21:23 -; author: Liam M. Healy -; modified: Sun Oct 8 2006 - 21:45 -;******************************************************** -;;; $Id: $ +;; Pareto distribution +;; Liam Healy, Sat Oct 8 2006 - 21:23 +;; Time-stamp: <2008-02-03 10:33:14EST pareto.lisp> +;; $Id: $ (in-package :gsl) @@ -13,43 +9,47 @@ "gsl_ran_pareto" (((generator generator) :pointer) (a :double) (b :double)) :c-return :double - :documentation - "A random variate from the Pareto distribution of order @var{a}. + :documentation ; FDL + "A random variate from the Pareto distribution of order a. The distribution function is p(x) dx = (a/b) / (x/b)^{a+1} dx - @math{x >= b}.") + x >= b.") (defun-gsl pareto-pdf (x a b) "gsl_ran_pareto_pdf" ((x :double) (a :double) (b :double)) :c-return :double - :documentation - "The probability density @math{p(x)} at @var{x} - for a Pareto distribution with exponent @var{a} and scale @var{b}, using + :documentation ; FDL + "The probability density p(x) at x + for a Pareto distribution with exponent a and scale b, using the formula given in #'pareto.") (defun-gsl pareto-P (x a b) "gsl_cdf_pareto_P" ((x :double) (a :double) (b :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{P(x)} for the Pareto distribution with exponent @var{a} and scale @var{b}.") + :documentation ; FDL + "The cumulative distribution functions + P(x) for the Pareto distribution with exponent a and scale b.") (defun-gsl pareto-Q (x a b) "gsl_cdf_pareto_Q" ((x :double) (a :double) (b :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{Q(x)} for the Pareto distribution with exponent @var{a} and scale @var{b}.") + :documentation ; FDL + "The cumulative distribution functions + Q(x) for the Pareto distribution with exponent a and scale b.") (defun-gsl pareto-Pinv (P a b) "gsl_cdf_pareto_Pinv" ((P :double) (a :double) (b :double)) :c-return :double - :documentation "The inverse cumulative distribution functions - @math{P(x)} for the Pareto distribution with exponent @var{a} and scale @var{b}.") + :documentation ; FDL + "The inverse cumulative distribution functions + P(x) for the Pareto distribution with exponent a and scale b.") (defun-gsl pareto-Qinv (Q a b) "gsl_cdf_pareto_Qinv" ((Q :double) (a :double) (b :double)) :c-return :double - :documentation "The inverse cumulative distribution functions - @math{Q(x)} for the Pareto distribution with exponent @var{a} and scale @var{b}.") + :documentation ; FDL + "The inverse cumulative distribution functions + Q(x) for the Pareto distribution with exponent a and scale b.") ;;; Examples and unit test (lisp-unit:define-test pareto @@ -59,11 +59,10 @@ "0.208882311615d+01" "0.268706924980d+01" "0.370340428797d+01" "0.270287443943d+01" "0.263177356639d+01") (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 collect - (pareto *rng-mt19937* 1.0d0 2.0d0))))) + (pareto rng 1.0d0 2.0d0))))) (lisp-unit:assert-first-fp-equal "0.511603440571d+00" (pareto-pdf 1.5d0 1.3d0 1.0d0)) diff --git a/random/poisson.lisp b/random/poisson.lisp index 70d5bf268d666732c39489e024e7ae5c69d2f40d..e12c9fa794cee1c63861c0a5067ae9b119141661 100644 --- a/random/poisson.lisp +++ b/random/poisson.lisp @@ -1,6 +1,6 @@ ;; Poisson distribution ;; Liam Healy, Sat Nov 25 2006 - 16:00 -;; Time-stamp: <2008-02-02 23:02:43EST poisson.lisp> +;; Time-stamp: <2008-02-03 09:50:16EST poisson.lisp> ;; $Id: $ (in-package :gsl) @@ -9,31 +9,33 @@ "gsl_ran_poisson" (((generator generator) :pointer) (mu :double)) :c-return :uint - :documentation - "A random integer from the Poisson distribution with mean @var{mu}. + :documentation ; FDL + "A random integer from the Poisson distribution with mean mu. The probability distribution for Poisson variates is p(k) = {\mu^k \over k!} \exp(-\mu) - @math{k >= 0}.") + k >= 0.") (defun-gsl poisson-pdf (k mu) "gsl_ran_poisson_pdf" ((k :uint) (mu :double)) :c-return :double - :documentation - "The probability @math{p(k)} of obtaining @var{k} - from a Poisson distribution with mean @var{mu} using the formula + :documentation ; FDL + "The probability p(k) of obtaining k + from a Poisson distribution with mean mu using the formula given in #'poisson.") (defun-gsl poisson-P (k mu) "gsl_cdf_poisson_P" ((k :uint) (mu :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{P(k)} for the Poisson distribution with parameter @var{mu}.") + :documentation ; FDL + "The cumulative distribution functions + P(k) for the Poisson distribution with parameter mu.") (defun-gsl poisson-Q (k mu) "gsl_cdf_poisson_Q" ((k :uint) (mu :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{Q(k)} for the Poisson distribution with parameter @var{mu}.") + :documentation ; FDL + "The cumulative distribution functions + Q(k) for the Poisson distribution with parameter mu.") ;;; Examples and unit test (lisp-unit:define-test poisson diff --git a/random/rayleigh-tail.lisp b/random/rayleigh-tail.lisp index a4671180fb9161d9ef803e100dae143de6eb11ce..2d6f00edac40fd54ee1ddfae654c1d83e17d3627 100644 --- a/random/rayleigh-tail.lisp +++ b/random/rayleigh-tail.lisp @@ -1,6 +1,6 @@ ;; Rayleigh tail distribution ;; Liam Healy, Sat Sep 30 2006 -;; Time-stamp: <2008-02-02 23:00:39EST rayleigh-tail.lisp> +;; Time-stamp: <2008-02-03 09:47:12EST rayleigh-tail.lisp> ;; $Id: $ (in-package :gsl) @@ -9,20 +9,20 @@ "gsl_ran_rayleigh_tail" (((generator generator) :pointer) (a :double) (sigma :double)) :c-return :double - :documentation + :documentation ; FDL "A random variate from the tail of the Rayleigh - distribution with scale parameter @var{sigma} and a lower limit of - @var{a}. The distribution is + distribution with scale parameter sigma and a lower limit of + a. The distribution is p(x) dx = {x \over \sigma^2} \exp ((a^2 - x^2) /(2 \sigma^2)) dx - for @math{x > a}.") + for x > a.") (defun-gsl rayleigh-tail-pdf (x a sigma) "gsl_ran_rayleigh_tail_pdf" ((x :double) (a :double) (sigma :double)) :c-return :double - :documentation - "The probability density @math{p(x)} at @var{x} - for a Rayleigh tail distribution with scale parameter @var{sigma} and - lower limit @var{a}, using the formula given in #'rayleigh-tail.") + :documentation ; FDL + "The probability density p(x) at x + for a Rayleigh tail distribution with scale parameter sigma and + lower limit a, using the formula given in #'rayleigh-tail.") ;;; Examples and unit test (lisp-unit:define-test rayleigh-tail diff --git a/random/rayleigh.lisp b/random/rayleigh.lisp index 693b7e8320cb7befe0f429be5d907187f0328ee2..7f66f1e3e1867cc9ca8e9e4889014c24e3ef83a5 100644 --- a/random/rayleigh.lisp +++ b/random/rayleigh.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: rayleigh.lisp -; description: Rayleigh distribution -; date: Sat Sep 30 2006 -; author: Liam M. Healy -; modified: Sat Sep 30 2006 - 19:41 -;******************************************************** -;;; $Id: $ +;; Rayleigh distribution +;; Liam Healy, Sat Sep 30 2006 +;; Time-stamp: <2008-02-03 11:29:48EST rayleigh.lisp> +;; $Id: $ (in-package :gsl) @@ -13,47 +9,51 @@ "gsl_ran_rayleigh" (((generator generator) :pointer) (sigma :double)) :c-return :double - :documentation + :documentation ; FDL "A random variate from the Rayleigh distribution with - scale parameter @var{sigma}. The distribution is + scale parameter sigma. The distribution is p(x) dx = {x \over \sigma^2} \exp(- x^2/(2 \sigma^2)) dx - for @math{x > 0}.") + for x > 0.") (defun-gsl rayleigh-pdf (x sigma) "gsl_ran_rayleigh_pdf" ((x :double) (sigma :double)) :c-return :double - :documentation - "The probability density @math{p(x)} at @var{x} - for a Rayleigh distribution with scale parameter @var{sigma}, using the + :documentation ; FDL + "The probability density p(x) at x + for a Rayleigh distribution with scale parameter sigma, using the formula given for #'rayleigh.") (defun-gsl rayleigh-P (x sigma) "gsl_cdf_rayleigh_P" ((x :double) (sigma :double)) :c-return :double - :documentation "The cumulative distribution function - @math{P(x)} for the Rayleigh distribution with scale - parameter @var{sigma}.") + :documentation ; FDL + "The cumulative distribution function + P(x) for the Rayleigh distribution with scale + parameter sigma.") (defun-gsl rayleigh-Q (x sigma) "gsl_cdf_rayleigh_Q" ((x :double) (sigma :double)) :c-return :double - :documentation "The cumulative distribution function - @math{Q(x)} for the Rayleigh distribution with scale - parameter @var{sigma}.") + :documentation ; FDL + "The cumulative distribution function + Q(x) for the Rayleigh distribution with scale + parameter sigma.") (defun-gsl rayleigh-Pinv (P sigma) "gsl_cdf_rayleigh_Pinv" ((P :double) (sigma :double)) :c-return :double - :documentation "The inverse cumulative distribution function - @math{P(x)} for the Rayleigh distribution with scale - parameter @var{sigma}.") + :documentation ; FDL + "The inverse cumulative distribution function + P(x)} for the Rayleigh distribution with scale + parameter sigma.") (defun-gsl rayleigh-Qinv (Q sigma) "gsl_cdf_rayleigh_Qinv" ((Q :double) (sigma :double)) :c-return :double - :documentation "The inverse cumulative distribution function - @math{Q(x)} for the Rayleigh distribution with scale - parameter @var{sigma}.") + :documentation ; FDL + "The inverse cumulative distribution function + Q(x) for the Rayleigh distribution with scale + parameter sigma.") ;;; Examples and unit test (lisp-unit:define-test rayleigh @@ -63,11 +63,10 @@ "0.294800354467d+01" "0.768510144246d+01" "0.111004981321d+02" "0.776103902005d+01" "0.740959915506d+01") (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 collect - (rayleigh *rng-mt19937* 10.0d0))))) + (rayleigh rng 10.0d0))))) (lisp-unit:assert-first-fp-equal "0.441248451292d+00" (rayleigh-pdf 0.5d0 1.0d0)) diff --git a/random/spherical-vector.lisp b/random/spherical-vector.lisp index 3b8147ebc80f808545a4b57128fe00b1ab802905..12a54f7ef9e308eb48f1e178a3a88cb157f3e3d7 100644 --- a/random/spherical-vector.lisp +++ b/random/spherical-vector.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: spherical-vector.lisp -; description: Spherical Vector distribution -; date: Sun Oct 22 2006 -; author: Liam M. Healy -; modified: Sun Oct 22 2006 - 21:58 -;******************************************************** -;;; $Id: $ +;; Spherical Vector distribution +;; Liam Healy, Sun Oct 22 2006 +;; Time-stamp: <2008-02-03 10:38:08EST spherical-vector.lisp> +;; $Id: $ (in-package :gsl) @@ -15,28 +11,28 @@ "gsl_ran_dir_2d" (((generator generator) :pointer) (x :double) (y :double)) :c-return :void - :documentation - "A random direction vector @math{v} = (@var{x},@var{y}) in + :documentation ; FDL + "A random direction vector v = (x,y) in two dimensions. The vector is normalized such that - @math{|v|^2 = x^2 + y^2 = 1}.") + |v|^2 = x^2 + y^2 = 1.") (defun-gsl direction-2d-trig-method (generator) "gsl_ran_dir_2d_trig_method" (((generator generator) :pointer) (x :double) (y :double)) :c-return :void - :documentation - "A random direction vector @math{v} = (@var{x},@var{y}) in + :documentation ; FDL + "A random direction vector v = (x,y) in two dimensions. The vector is normalized such that - @math{|v|^2 = x^2 + y^2 = 1}. Uses trigonometric functions.") + |v|^2 = x^2 + y^2 = 1. Uses trigonometric functions.") (defun-gsl direction-3d (generator) "gsl_ran_dir_3d" (((generator generator) :pointer) (x :double) (y :double) (z :double)) :c-return :void - :documentation - "A random direction vector @math{v} = - (@var{x},@var{y},@var{z}) in three dimensions. The vector is normalized - such that @math{|v|^2 = x^2 + y^2 + z^2 = 1}. The method employed is + :documentation ; FDL + "A random direction vector v = + (x,y,z) in three dimensions. The vector is normalized + such that |v|^2 = x^2 + y^2 + z^2 = 1. The method employed is due to Robert E. Knop (CACM 13, 326 (1970)), and explained in Knuth, v2, 3rd ed, p136. It uses the surprising fact that the distribution projected along any axis is actually uniform (this is only true for 3 @@ -47,10 +43,10 @@ (((generator generator) :pointer) ((dim0 x) :size) ((gsl-array x) :pointer)) :c-return :void :return (x) - :documentation - "A random direction vector @math{v = (x_1,x_2,...,x_n)} in @var{n} dimensions, + :documentation ; FDL + "A random direction vector v = (x_1,x_2,...,x_n) in n dimensions, where n is the length of the vector x passed in. The vector is normalized such that - @math{|v|^2 = x_1^2 + x_2^2 + ... + x_n^2 = 1}. The method + |v|^2 = x_1^2 + x_2^2 + ... + x_n^2 = 1. The method uses the fact that a multivariate gaussian distribution is spherically symmetric. Each component is generated to have a gaussian distribution, and then the components are normalized. The method is described by @@ -65,29 +61,26 @@ "0.457266229462d+00" "0.889329857473d+00" "-0.463256161598d+00" "-0.886224423462d+00") (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 4 append - (multiple-value-list (direction-2d *rng-mt19937*)))))) + (multiple-value-list (direction-2d rng)))))) (lisp-unit:assert-equal '("0.999998683521d+00" "-0.162263876311d-02" "0.520301010608d+00" "0.853982937980d+00" "-0.203512053104d+00" "0.979072440753d+00" "0.945475322749d+00" "-0.325693742761d+00" "0.115000339166d+00" "0.993365452385d+00") (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 4 append - (multiple-value-list (direction-2d-trig-method *rng-mt19937*)))))) + (multiple-value-list (direction-2d-trig-method rng)))))) (lisp-unit:assert-equal '("-0.912992575045d-01" "0.187821853572d+00" "0.977950610665d+00" "-0.905118296156d+00" "-0.506837644858d-01" "-0.422127973465d+00" "0.139937665360d+00" "0.838546262052d+00" "-0.526552576873d+00") (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 2 append - (multiple-value-list (direction-3d *rng-mt19937*))))))) + (multiple-value-list (direction-3d rng))))))) diff --git a/random/tdist.lisp b/random/tdist.lisp index 667f8f7c89c0e7ae442a998c780e2759371e13b1..9a88f88ca7b2357b83c6f51439ce399a163c5220 100644 --- a/random/tdist.lisp +++ b/random/tdist.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: tdist.lisp -; description: tdist distribution -; date: Sat Oct 7 2006 - 16:13 -; author: Liam M. Healy -; modified: Sun Oct 8 2006 - 16:44 -;******************************************************** -;;; $Id: $ +;; Tdist distribution +;; Liam Healy, Sat Oct 7 2006 - 16:13 +;; Time-stamp: <2008-02-03 11:27:51EST tdist.lisp> +;; $Id: $ (in-package :gsl) @@ -13,44 +9,48 @@ "gsl_ran_tdist" (((generator generator) :pointer) (nu :double)) :c-return :double - :documentation + :documentation ; FDL "A random variate from the Student t-distribution. The distribution function is, p(x) dx = {\Gamma((\nu + 1)/2) \over \sqrt{\pi \nu} \Gamma(\nu/2)} (1 + x^2/\nu)^{-(\nu + 1)/2} dx - for @math{-\infty < x < +\infty}.") + for -\infty < x < +\infty.") (defun-gsl tdist-pdf (x nu) "gsl_ran_tdist_pdf" ((x :double) (nu :double)) :c-return :double - :documentation - "The probability density @math{p(x)} at @var{x} - for a t-distribution with @var{nu} degrees of freedom, using the formula + :documentation ; FDL + "The probability density p(x) at x + for a t-distribution with nu degrees of freedom, using the formula given in #'tdist.") (defun-gsl tdist-P (x nu) "gsl_cdf_tdist_P" ((x :double) (nu :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{P(x)} for the tdist distribution with @var{nu} degrees of freedom.") + :documentation ; FDL + "The cumulative distribution functions + P(x) for the tdist distribution with nu degrees of freedom.") (defun-gsl tdist-Q (x nu) "gsl_cdf_tdist_Q" ((x :double) (nu :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{Q(x)} for the tdist distribution with @var{nu} degrees of freedom.") + :documentation ; FDL + "The cumulative distribution functions + Q(x) for the tdist distribution with nu degrees of freedom.") (defun-gsl tdist-Pinv (P nu) "gsl_cdf_tdist_Pinv" ((P :double) (nu :double)) :c-return :double - :documentation "The inverse cumulative distribution functions - @math{P(x)} for the tdist distribution with @var{nu} degrees of freedom.") + :documentation ; FDL + "The inverse cumulative distribution functions + P(x) for the tdist distribution with nu degrees of freedom.") (defun-gsl tdist-Qinv (Q nu) "gsl_cdf_tdist_Qinv" ((Q :double) (nu :double)) :c-return :double - :documentation "The inverse cumulative distribution functions - @math{Q(x)} for the tdist distribution with @var{nu} degrees of freedom.") + :documentation ; FDL + "The inverse cumulative distribution functions + Q(x) for the tdist distribution with nu degrees of freedom.") ;;; Examples and unit test (lisp-unit:define-test tdist @@ -59,12 +59,11 @@ "-0.160088628258d+01" "-0.170109355058d+01" "-0.437095974981d-01" "0.127611592766d+00" "-0.197312182555d-01" "-0.653466611720d+00" "0.203577132452d+00" "0.177650300478d+01") - (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) - (loop for i from 0 to 10 - collect - (tdist *rng-mt19937* 10.0d0))))) + (letm ((rng (random-number-generator *mt19937* 0))) + (rng-set rng 0) + (loop for i from 0 to 10 + collect + (tdist rng 10.0d0))))) (lisp-unit:assert-first-fp-equal "0.254647908947d+00" (tdist-pdf 0.5d0 1.0d0)) diff --git a/random/weibull.lisp b/random/weibull.lisp index 44da1a63bf64bdbb804ec7776e945ff909f37ac7..c857e81907b5f7abfa642977c29c2eae6841b1f4 100644 --- a/random/weibull.lisp +++ b/random/weibull.lisp @@ -1,11 +1,7 @@ -;******************************************************** -; file: weibull.lisp -; description: Weibull distribution -; date: Sun Oct 22 2006 -; author: Liam M. Healy -; modified: Sat Nov 11 2006 - 20:59 -;******************************************************** -;;; $Id: $ +;; Weibull distribution +;; Liam Healy, Sun Oct 22 2006 +;; Time-stamp: <2008-02-03 11:15:10EST weibull.lisp> +;; $Id: $ (in-package :gsl) @@ -13,42 +9,46 @@ "gsl_ran_weibull" (((generator generator) :pointer) (a :double) (b :double)) :c-return :double - :documentation + :documentation ; FDL "A random variate from the Weibull distribution. The distribution function is p(x) dx = {b \over a^b} x^{b-1} \exp(-(x/a)^b) dx - for @math{x >= 0}.") + for x >= 0.") (defun-gsl weibull-pdf (x a b) "gsl_ran_weibull_pdf" ((x :double) (a :double) (b :double)) :c-return :double - :documentation - "The probability density @math{p(x)} at @var{x} - for a Weibull distribution with scale @var{a} and exponent @var{b}, + :documentation ; FDL + "The probability density p(x) at x + for a Weibull distribution with scale a and exponent b, using the formula given in #'weibull.") (defun-gsl weibull-P (x a b) "gsl_cdf_weibull_P" ((x :double) (a :double) (b :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{P(x)} for the Weibull distribution with scale @var{a} and exponent @var{b}.") + :documentation ; FDL + "The cumulative distribution functions + P(x) for the Weibull distribution with scale a and exponent b.") (defun-gsl weibull-Q (x a b) "gsl_cdf_weibull_Q" ((x :double) (a :double) (b :double)) :c-return :double - :documentation "The cumulative distribution functions - @math{Q(x)} for the Weibull distribution with scale @var{a} and exponent @var{b}.") + :documentation ; FDL + "The cumulative distribution functions + Q(x) for the Weibull distribution with scale a and exponent b.") (defun-gsl weibull-Pinv (P a b) "gsl_cdf_weibull_Pinv" ((P :double) (a :double) (b :double)) :c-return :double - :documentation "The inverse cumulative distribution functions - @math{P(x)} for the Weibull distribution scale @var{a} and exponent @var{b}.") + :documentation ; FDL + "The inverse cumulative distribution functions + P(x) for the Weibull distribution scale a and exponent b.") (defun-gsl weibull-Qinv (Q a b) "gsl_cdf_weibull_Qinv" ((Q :double) (a :double) (b :double)) :c-return :double - :documentation "The inverse cumulative distribution functions - @math{Q(x)} for the Weibull distribution exponent @var{a} and scale @var{b}.") + :documentation ; FDL + "The inverse cumulative distribution functions + Q(x) for the Weibull distribution exponent a and scale b.") ;;; Examples and unit test (lisp-unit:define-test weibull @@ -58,11 +58,10 @@ "0.208455329740d+00" "0.543418734407d+00" "0.784923750377d+00" "0.548788332013d+00" "0.523937780842d+00") (lisp-unit:fp-sequence - (progn - (rng-set *rng-mt19937* 0) + (letm ((rng (random-number-generator *mt19937* 0))) (loop for i from 0 to 10 collect - (weibull *rng-mt19937* 1.0d0 2.0d0))))) + (weibull rng 1.0d0 2.0d0))))) (lisp-unit:assert-first-fp-equal "0.242631749722d+00" (weibull-pdf 1.5d0 1.3d0 1.0d0))