Raymond Toy pushed to branch rtoy-xoro-default at cmucl / cmucl

Commits:

2 changed files:

Changes:

  • src/code/rand-xoroshiro.lisp
    ... ... @@ -8,7 +8,8 @@
    8 8 
     ;;; **********************************************************************
    9 9 
     ;;;
    10 10 
     ;;; Support for the xoroshiro128+ random number generator by David
    11 
    -;;; Blackman and Sebastiano Vigna (vigna@acm.org)
    11 
    +;;; Blackman and Sebastiano Vigna (vigna@acm.org). See
    12 
    +;;; http://xoroshiro.di.unimi.it/.
    12 13
    13 14 
     (in-package "LISP")
    14 15 
     (intl:textdomain "cmucl")
    ... ... @@ -47,6 +48,18 @@
    47 48 
       (let ((state (or state (make-array 2 :element-type 'double-float)))
    48 49 
     	(splitmix-state (ldb (byte 64 0) seed)))
    49 50 
         (flet ((splitmix64 ()
    51 
    +	     ;; See http://xoroshiro.di.unimi.it/splitmix64.c for the
    52 
    +	     ;; definitive reference.  The basic algorithm, where x is
    53 
    +	     ;; the 64-bit state of the generator,:
    54 
    +	     ;;
    55 
    +	     ;;   uint64_t z = (x += 0x9e3779b97f4a7c15);
    56 
    +	     ;;   z = (z ^ (z >> 30)) * 0xbf58476d1ce4e5b9;
    57 
    +	     ;;   z = (z ^ (z >> 27)) * 0x94d049bb133111eb;
    58 
    +	     ;;   return z ^ (z >> 31);
    59 
    +	     ;;
    60 
    +	     ;; This is only used occasionally for initializing the
    61 
    +	     ;; RNG, so this is a very straight-forward
    62 
    +	     ;; implementation.
    50 63 
     	     (let ((z (setf splitmix-state
    51 64 
     			    (ldb (byte 64 0) (+ splitmix-state #x9e3779b97f4a7c15)))))
    52 65 
     	       (declare (type (unsigned-byte 64) z))
    ... ... @@ -192,8 +205,8 @@
    192 205
    193 206 
     ;;;; Random entries:
    194 207
    195 
    -;;#+x86
    196 
    -;;(declaim (inline xoroshiro-next))
    208 
    +#+x86
    209 
    +(declaim (inline xoroshiro-gen))
    197 210 
     #+x86
    198 211 
     (defun xoroshiro-gen (state)
    199 212 
       (declare (type (simple-array double-float (2)) state)
    ... ... @@ -204,7 +217,31 @@
    204 217 
     (defun xoroshiro-gen (state)
    205 218 
       (declare (type (simple-array double-float (2)) state)
    206 219 
     	   (optimize (speed 3) (safety 0)))
    220 
    +  ;; Portable implemenation of the xoroshiro128+ generator. See
    221 
    +  ;; http://xoroshiro.di.unimi.it/xoroshiro128plus.c for the
    222 
    +  ;; definitive definition.
    223 
    +  ;;
    224 
    +  ;; uint64_t s[2];
    225 
    +  ;;
    226 
    +  ;; static inline uint64_t rotl(const uint64_t x, int k) {
    227 
    +  ;; 	return (x << k) | (x >> (64 - k));
    228 
    +  ;; }
    229 
    +  ;;
    230 
    +  ;; uint64_t next(void) {
    231 
    +  ;; 	const uint64_t s0 = s[0];
    232 
    +  ;; 	uint64_t s1 = s[1];
    233 
    +  ;; 	const uint64_t result = s0 + s1;
    234 
    +  ;;
    235 
    +  ;; 	s1 ^= s0;
    236 
    +  ;; 	s[0] = rotl(s0, 55) ^ s1 ^ (s1 << 14); // a, b
    237 
    +  ;; 	s[1] = rotl(s1, 36); // c
    238 
    +  ;;
    239 
    +  ;; 	return result;
    240 
    +  ;; }
    241 
    +  ;;
    207 242 
       (flet ((rotl-55 (x1 x0)
    243 
    +	   ;; Rotate [x1|x0] left 55 bits, returning the result as two
    244 
    +	   ;; values.
    208 245 
     	   (declare (type (unsigned-byte 32) x0 x1)
    209 246 
     		    (optimize (speed 3) (safety 0)))
    210 247 
     	   ;; x << 55
    ... ... @@ -218,6 +255,8 @@
    218 255 
     	       (values (logior sl55-h sr9-h)
    219 256 
     		       (logior sl55-l sr9-l)))))
    220 257 
     	 (rotl-36 (x1 x0)
    258 
    +	   ;; Rotate [x1|x0] left 36 bits, returning the result as two
    259 
    +	   ;; values.
    221 260 
     	   (declare (type (unsigned-byte 32) x0 x1)
    222 261 
     		    (optimize (speed 3) (safety 0)))
    223 262 
     	   ;; x << 36
    ... ... @@ -230,6 +269,8 @@
    230 269 
     	       (values (logior sl36-h sr28-h)
    231 270 
     		       sr28-l))))
    232 271 
     	 (shl-14 (x1 x0)
    272 
    +	   ;; Shift [x1|x0] left by 14 bits, returning the result as
    273 
    +	   ;; two values.
    233 274 
     	   (declare (type (unsigned-byte 32) x1 x0)
    234 275 
     		    (optimize (speed 3) (safety 0)))
    235 276 
     	   (values (ldb (byte 32 0)
    ... ... @@ -248,6 +289,9 @@
    248 289 
     	  (s1-1 0)
    249 290 
     	  (s1-0 0))
    250 291 
           (declare (type (unsigned-byte 32) s0-1 s0-0 s1-1 s1-0))
    292 
    +      ;; Load the state to s0 and s1. s0-1 is the high 32-bit part and
    293 
    +      ;; s0-0 is the low 32-bit part of the 64-bit value.  Similarly
    294 
    +      ;; for s1.
    251 295 
           (multiple-value-bind (x1 x0)
    252 296 
     	  (kernel:double-float-bits (aref state 0))
    253 297 
     	(setf s0-1 (ldb (byte 32 0) x1)
    ... ... @@ -257,6 +301,7 @@
    257 301 
     	(setf s1-1 (ldb (byte 32 0) x1)
    258 302 
     	      s1-0 x0))
    259 303
    304 
    +      ;; Compute the 64-bit random value: s0 + s1
    260 305 
           (multiple-value-prog1
    261 306 
     	  (multiple-value-bind (sum-0 c)
    262 307 
     	      (bignum::%add-with-carry s0-0 s1-0 0)

  • src/compiler/x86/arith.lisp
    ... ... @@ -1835,51 +1835,8 @@
    1835 1835 
     	   (give-up)))))
    1836 1836
    1837 1837 
     (in-package "VM")
    1838 
    -#+nil
    1839 
    -(progn
    1840 
    -(defknown xoroshiro-next (double-float double-float)
    1841 
    -  (values (unsigned-byte 32) (unsigned-byte 32) double-float double-float)
    1842 
    -  (movable))
    1843 
    -
    1844 
    -(define-vop (xoroshiro-next)
    1845 
    -  (:policy :fast-safe)
    1846 
    -  (:translate xoroshiro-next)
    1847 
    -  (:args (old-s1 :scs (double-reg) :to (:result 3))
    1848 
    -	 (old-s0 :scs (double-reg) :to (:result 3)))
    1849 
    -  (:arg-types double-float double-float)
    1850 
    -  (:results (r1 :scs (unsigned-reg))
    1851 
    -	    (r0 :scs (unsigned-reg))
    1852 
    -	    (s1 :scs (double-reg))
    1853 
    -	    (s0 :scs (double-reg)))
    1854 
    -  (:result-types unsigned-num unsigned-num double-float double-float)
    1855 
    -  (:temporary (:sc double-reg) t0)
    1856 
    -  (:generator 10
    1857 
    -    (inst movapd t0 old-s0)
    1858 
    -    (inst paddq t0 old-s1)		; t0 = old-s0 + old-s1
    1859 
    -    (inst movd r0 t0)			; r0 = low 32-bits of t0
    1860 
    -    (inst psrlq t0 32)
    1861 
    -    (inst movd r1 t0)			; r1 = high 32-bits of t0
    1862 
    -    ;; s1 ^= s0
    1863 
    -    (inst movapd s1 old-s1)		; s1 = old-s1
    1864 
    -    (inst xorpd s1 old-s0)		; s1 = old-s0 ^ old-s1
    1865 
    -    ;; rotl(s0, 55) = s0 << 55 | (s0 >> 9)
    1866 
    -    (inst movapd s0 old-s0)		; s0 = old-s0
    1867 
    -    (inst movapd t0 old-s0)		; t0 = old-s0
    1868 
    -    (inst psllq s0 55)			; s0 = s0 << 55
    1869 
    -    (inst psrlq t0 9)			; t0 = s0 >> 9
    1870 
    -    (inst orpd s0 t0)	                ; s0 = rotl(s0,55) = s0 << 55 | s0 >> 9
    1871 
    -    (inst xorpd s0 s1)			; s0 = rotl(s0,55) ^ s1
    1872 
    -    (inst movapd t0 s1)			; t0 = s1
    1873 
    -    (inst psllq t0 14)			; t0 = s1 << 14
    1874 
    -    (inst xorpd s0 t0)		        ; s0 = rotl(s0,55) ^ s1 ^ (s1 << 14)
    1875 
    -    (inst movapd t0 s1)			; t0 = s1
    1876 
    -    (inst psllq t0 36)			; t0 = s1 << 36
    1877 
    -    (inst psrlq s1 28)			; s1 = s1 >> 28
    1878 
    -    (inst orpd s1 t0)			; s1 = rotl(new-s1, 36)
    1879 
    -
    1880 
    -    ))
    1881 
    -)
    1882 1838
    1839 
    +#+random-xoroshiro
    1883 1840 
     (progn
    1884 1841 
     (defknown xoroshiro-next ((simple-array double-float (2)))
    1885 1842 
       (values (unsigned-byte 32) (unsigned-byte 32))

    • —
    View it on GitLab.
    You're receiving this email because of your account on gitlab.common-lisp.net. If you'd like to receive fewer emails, you can adjust your notification settings.