## Shellsort With Three Increments

### December 15, 2015

We begin by implementing shellsort, giving the list of gaps as a parameter; for simplicity, we assume we are sorting numbers:

(define (shellsort vec gaps) (let ((n (vector-length vec))) (do ((gaps gaps (cdr gaps))) ((null? gaps) vec) (let ((g (car gaps))) (do ((i g (+ i 1))) ((<= n i)) (let ((t (vector-ref vec i))) (do ((j i (- j g))) ((or (< j g) (<= (vector-ref vec (- j g)) t)) (vector-set! vec j t)) (vector-set! vec j (vector-ref vec (- j g))))))))))

Next we implement the selection of the gap sequence; *h* is just the square root of *n*, *g* starts at the square root of *h*, then decrements until it is coprime to *h*:

(define (gaps n) (define (coprime? x y) (= (gcd x y) 1)) (let* ((h (isqrt n)) (g (isqrt h))) (let loop ((g g)) (if (coprime? h g) (list h g 1) (loop (- g 1))))))

For our timing tests we need a way to create vectors of *n* random numbers:

(define (rand-vector n) (let ((vec (make-vector n))) (do ((i 0 (+ i 1))) ((= i n) vec) (vector-set! vec i (randint 1000000)))))

We’ll run the timing tests *k* times on random vectors of length *n*, reporting the average cpu time for each of the *k* sorts, in milliseconds:

(define (time-sort n k) (let loop ((m k) (sum 0)) (if (zero? m) (round (/ sum k)) (let ((vec (rand-vector n))) (let ((start (cpu-time))) (let ((ignore (shellsort vec (gaps n)))) (loop (- m 1) (+ sum (- (cpu-time) start)))))))))

Now we are ready for some timings:

> (time-sort 10000 10) 90 > (time-sort 20000 10) 253 > (time-sort 40000 10) 690 > (time-sort 80000 10) 1995 > (time-sort 160000 6) 5530 > (time-sort 320000 2) 15592

As we double the size of the input, the timings increase as a factor of 2.81, 2.73, 2.89, 2.77 and 2.82, which is clearly somewhere between linear and quadratic. Janson and Knuth’s conjecture of O(*n*^{3/2}) leads to a doubling factor of √8 ≈ 2.82, so we will declare the conjecture correct (I, for one, never doubted it).

We used isqrt and the random number generator from the Standard Prelude. You can run the program at http://ideone.com/LMEqyf.

In Python. For N > 1000 the conjecture seems to hold.