fork download
  1. ; shellsort with three increments
  2.  
  3. (define (isqrt n)
  4. (if (not (and (positive? n) (integer? n)))
  5. (error 'isqrt "must be positive integer")
  6. (let loop ((x n))
  7. (let ((y (quotient (+ x (quotient n x)) 2)))
  8. (if (< y x) (loop y) x)))))
  9.  
  10. (define rand #f)
  11. (define randint #f)
  12. (let ((two31 #x80000000) (a (make-vector 56 -1)) (fptr #f))
  13. (define (mod-diff x y) (modulo (- x y) two31)) ; generic version
  14. ; (define (mod-diff x y) (logand (- x y) #x7FFFFFFF)) ; fast version
  15. (define (flip-cycle)
  16. (do ((ii 1 (+ ii 1)) (jj 32 (+ jj 1))) ((< 55 jj))
  17. (vector-set! a ii (mod-diff (vector-ref a ii) (vector-ref a jj))))
  18. (do ((ii 25 (+ ii 1)) (jj 1 (+ jj 1))) ((< 55 ii))
  19. (vector-set! a ii (mod-diff (vector-ref a ii) (vector-ref a jj))))
  20. (set! fptr 54) (vector-ref a 55))
  21. (define (init-rand seed)
  22. (let* ((seed (mod-diff seed 0)) (prev seed) (next 1))
  23. (vector-set! a 55 prev)
  24. (do ((i 21 (modulo (+ i 21) 55))) ((zero? i))
  25. (vector-set! a i next) (set! next (mod-diff prev next))
  26. (set! seed (+ (quotient seed 2) (if (odd? seed) #x40000000 0)))
  27. (set! next (mod-diff next seed)) (set! prev (vector-ref a i)))
  28. (flip-cycle) (flip-cycle) (flip-cycle) (flip-cycle) (flip-cycle)))
  29. (define (next-rand)
  30. (if (negative? (vector-ref a fptr)) (flip-cycle)
  31. (let ((next (vector-ref a fptr))) (set! fptr (- fptr 1)) next)))
  32. (define (unif-rand m)
  33. (let ((t (- two31 (modulo two31 m))))
  34. (let loop ((r (next-rand)))
  35. (if (<= t r) (loop (next-rand)) (modulo r m)))))
  36. (init-rand 19380110) ; happy birthday donald e knuth
  37. (set! rand (lambda seed
  38. (cond ((null? seed) (/ (next-rand) two31))
  39. ((eq? (car seed) 'get) (cons fptr (vector->list a)))
  40. ((eq? (car seed) 'set) (set! fptr (caadr seed))
  41. (set! a (list->vector (cdadr seed))))
  42. (else (/ (init-rand (modulo (numerator
  43. (inexact->exact (car seed))) two31)) two31)))))
  44. (set! randint (lambda args
  45. (cond ((null? (cdr args))
  46. (if (< (car args) two31) (unif-rand (car args))
  47. (floor (* (next-rand) (car args)))))
  48. ((< (car args) (cadr args))
  49. (let ((span (- (cadr args) (car args))))
  50. (+ (car args)
  51. (if (< span two31) (unif-rand span)
  52. (floor (* (next-rand) span))))))
  53. (else (let ((span (- (car args) (cadr args))))
  54. (- (car args)
  55. (if (< span two31) (unif-rand span)
  56. (floor (* (next-rand) span))))))))))
  57.  
  58. (define (shellsort lt? vec gaps)
  59. (let ((n (vector-length vec)))
  60. (do ((gaps gaps (cdr gaps))) ((null? gaps) vec)
  61. (let ((g (car gaps)))
  62. (do ((i g (+ i 1))) ((<= n i))
  63. (let ((t (vector-ref vec i)))
  64. (do ((j i (- j g)))
  65. ((or (< j g)
  66. (<= (vector-ref vec (- j g)) t))
  67. (vector-set! vec j t))
  68. (vector-set! vec j
  69. (vector-ref vec (- j g))))))))))
  70.  
  71. (define (gaps n)
  72. (define (coprime? x y) (= (gcd x y) 1))
  73. (let* ((h (isqrt n)) (g (isqrt h)))
  74. (let loop ((g g))
  75. (if (coprime? h g)
  76. (list h g 1)
  77. (loop (- g 1))))))
  78.  
  79. (define (rand-vector n)
  80. (let ((vec (make-vector n)))
  81. (do ((i 0 (+ i 1))) ((= i n) vec)
  82. (vector-set! vec i (randint 1000000)))))
  83.  
  84. (display (shellsort < (rand-vector 50) (gaps 50)))
Success #stdin #stdout 0.09s 8968KB
stdin
Standard input is empty
stdout
#(22908 32377 38095 48883 55898 61116 67109 96346 111793 117753 223602 226214 229380 230438 235029 247050 275421 301018 304574 311250 312471 320438 388790 399788 448323 450160 506777 551994 570609 580564 612474 632064 632264 632947 681835 712382 725693 767224 781530 815609 862383 870840 874999 885689 894419 901955 940939 958342 970538 970699)