## 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(n3/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.

Pages: 1 2

### One Response to “Shellsort With Three Increments”

1. Paul said

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

```from fractions import gcd
from math import sqrt
from timeit import repeat

def shellsort(array, gaps=None):
n = len(array)
if gaps is None:
h, g = round(sqrt(n)), round(sqrt(sqrt(n)))
while gcd(h, g) != 1:
h += 1
gaps = (h, g, 1)
for gap in gaps:
for i in range(gap, n):
val = array[i]
j = i
while j >= gap and array[j - gap] > val:
array[j] = array[j - gap]
j -= gap
array[j] = val

if __name__ == "__main__":
for n in range(7):
N = 10 ** n
setup = "import random;from __main__ import shellsort;"
setup += "data = cdata = list(range(" + str(N) + "));"
setup += "random.shuffle(data)"
stmt = "shellsort(data)"
print("{:8d} {:12g} {:12g}".format(
N, N ** 1.5 / 3000000, min(repeat(stmt, setup, number=6))))
"""     N     cN^3/2      time
1  3.33333e-07  3.07887e-05
10  1.05409e-05  0.000104681
100  0.000333333  0.000972409
1000    0.0105409    0.0162739
10000     0.333333     0.364905
100000      10.5409      10.0037
1000000      333.333      322.088
"""
```