1 /** Fixed Length Sorting via Sorting Networks. 2 3 See_Also: http://forum.dlang.org/post/ne5m62$1gu5$1@digitalmars.com 4 See_Also: http://cpansearch.perl.org/src/JGAMBLE/Algorithm-Networksort-1.30/lib/Algorithm/Networksort.pm 5 See_Also: http://www.angelfire.com/blog/ronz/Articles/999SortingNetworksReferen.html 6 7 License: $(WEB boost.org/LICENSE_1_0.txt, Boost License 1.0). 8 9 TODO see some sizes are not supported, we should not have holes. 10 Use http://www.angelfire.com/blog/ronz/Articles/999SortingNetworksReferen.html 11 12 TODO Sometimes the sort routine gets too bulky. Suggestion: also define 13 networks for `medianOfUpTo` and `medianExactly`, then use them in a 14 quicksort manner - first compute the median to segregate values in 15 below/over the median, then make two calls to `sortExactly!(n / 2)`. That 16 way you get to sort n values with median of n values (smaller and simpler) 17 and two sorts of n / 2 values. 18 19 TODO Stability of equal elements: Need template parameter `equalityStability`? Scalar builtin values are always stable. 20 21 TODO There should be a notion of at what point the networks become too bulky 22 to be fast - 6-16 may be the limit. 23 24 TODO There is a nice peephole optimization you could make. Consider: 25 26 r.conditionalSwap!("a < b", less, 0, 1, 1, 2)(r); 27 28 which needs to do 29 30 if (r[1] < r[0]) r.swapAt(0, 1); 31 if (r[2] < r[1]) r.swapAt(1, 2); 32 33 For types with no elaborate copy/assignment, it's more efficient to use a 34 "hole"-based approach - assign the first element to a temporary and then 35 consider it a hole that you fill, then leaving another hole: 36 37 if (r[1] < r[0]) 38 if (r[2] < r[0]) r.swapAt(0, 1, 2); 39 else r.swapAt(0, 1); 40 else 41 if (r[2] < r[1]) r.swapAt(1, 2); 42 43 with swapAt with three argument having this definition: 44 45 auto t = r[a]; r[a] = r[b]; r[b] = r[c]; r[c] = t; 46 47 i.e. create a temporary (which creates a "hole" in the array) then fill it 48 leaving another hole etc., until the last hole is filled with the temporary. 49 */ 50 module nxt.sortn; 51 52 import std.meta : allSatisfy; 53 import std.traits : isIntegral; 54 import std.range.primitives : isRandomAccessRange; 55 56 version(unittest) import std.algorithm.comparison : equal; 57 58 /** Static Iota. 59 TODO Move to Phobos std.range. 60 */ 61 template iota(size_t from, size_t to) 62 if (from <= to) 63 { 64 alias iota = iotaImpl!(to - 1, from); 65 } 66 private template iotaImpl(size_t to, size_t now) 67 { 68 import std.meta : AliasSeq; 69 static if (now >= to) { alias iotaImpl = AliasSeq!(now); } 70 else { alias iotaImpl = AliasSeq!(now, iotaImpl!(to, now + 1)); } 71 } 72 73 /** Conditionally pairwise sort elements of `Range` `r` at `indexes` using 74 comparison predicate `less`. 75 76 TODO Perhaps defines as 77 78 template conditionalSwap(indexes...) 79 { 80 void conditionalSwap(less = "a < b", Range)(Range r) 81 { 82 } 83 } 84 85 instead. 86 */ 87 void conditionalSwap(alias less = "a < b", Range, indexes...)(Range r) 88 if (isRandomAccessRange!Range && 89 allSatisfy!(isIntegral, typeof(indexes)) && 90 indexes.length && 91 (indexes.length & 1) == 0) // even number of indexes 92 { 93 import std.algorithm.mutation : swapAt; 94 import std.functional : binaryFun; 95 import nxt.static_iota : iota; 96 foreach (const i; iota!(0, indexes.length / 2)) 97 { 98 const j = indexes[2*i]; 99 const k = indexes[2*i + 1]; 100 101 static assert(j >= 0, "First part of index pair " ~ i.stringof ~ " is negative"); 102 static assert(k >= 0, "Second part of index pair " ~ i.stringof ~ " is negative"); 103 104 if (!binaryFun!less(r[j], r[k])) 105 { 106 r.swapAt(j, k); 107 } 108 } 109 } 110 111 /** Largest length supported by network sort `networkSortUpTo`. */ 112 enum networkSortMaxLength = 22; 113 114 /** Sort at most the first `n` elements of `r` using comparison `less` using 115 * a networking sort. 116 * 117 * Note: Sorting networks are not unique, not even optimal solutions. 118 * 119 * See_Also: http://stackoverflow.com/questions/3903086/standard-sorting-networks-for-small-values-of-n 120 * See_Also: http://www.cs.brandeis.edu/~hugues/sorting_networks.html 121 */ 122 auto networkSortUpTo(uint n, alias less = "a < b", Range)(Range r) 123 if (n >= 2 && 124 isRandomAccessRange!Range) 125 in 126 { 127 assert(r.length >= n); 128 } 129 do 130 { 131 auto s = r[0 .. n]; 132 133 static if (n < 2) 134 { 135 return; // no sorting needed 136 } 137 138 static if (n == 2) 139 { 140 s.conditionalSwap!(less, Range, 0, 1); 141 } 142 else static if (n == 3) 143 { 144 s.conditionalSwap!(less, Range, 145 0,1, 146 1,2, 147 0,1); 148 } 149 else static if (n == 4) 150 { 151 s.conditionalSwap!(less, Range, 152 0,1, 2,3, // 2 in parallel 153 0,2, 1,3, // 2 in parallel 154 1,2); 155 } 156 else static if (n == 5) 157 { 158 s.conditionalSwap!(less, Range, 159 0,1, 3,4, // 2 in parallel 160 0,2, 161 0,3, 1,2, // 2 in parallel 162 2,3, 1,4, // 2 in parallel 163 1,2, 3,4); // 2 in parallel 164 } 165 else static if (n == 6) 166 { 167 s.conditionalSwap!(less, Range, 168 0,1, 2,3, 4,5, // 3 in parallel 169 0,2, 1,4, 3,5, // 3 in parallel 170 0,1, 2,3, 4,5, // 3 in parallel 171 1,2, 3,4, // 2 in parallel 172 2,3); 173 } 174 else static if (n == 8) // Bitonic Sorter. 6-steps: https://en.wikipedia.org/wiki/Bitonic_sorter 175 { 176 s.conditionalSwap!(less, Range, 177 0,1, 2,3, 4,5, 6,7, 178 0,3, 4,7, 1,2, 5,6, 179 0,1, 2,3, 4,5, 6,7, 180 0,7, 1,6, 2,5, 3,4, 181 0,2, 1,3, 4,6, 5,7, 182 0,1, 2,3, 4,5, 6,7); 183 } 184 else static if (n == 9) // R. W. Floyd. 185 { 186 s.conditionalSwap!(less, Range, 187 0,1, 3,4, 6,7, 188 1,2, 4,5, 7,8, 189 0,1, 3,4, 6,7, 190 0,3, 191 3,6, 192 0,3, 1,4, 193 4,7, 194 1,4, 2,5, 195 5,8, 196 1,3, 2,5, 197 2,6, 5,7, 198 4,6, 199 2,4, 200 2,3, 5,6); 201 } 202 else static if (n == 10) // A. Waksman. 203 { 204 s.conditionalSwap!(less, Range, 205 0,5, 1,6, 2,7, 3,8, 4,9, 206 0,3, 1,4, 5,8, 6,9, 207 0,2, 3,6, 7,9, 208 0,1, 2,4, 5,7, 8,9, 209 1,2, 3,5, 4,6, 7,8, 210 1,3, 2,5, 4,7, 6,8, 211 2,3, 6,7, 212 3,4, 5,6, 213 4,5); 214 } 215 else static if (n == 11) // 12-input by Shapiro and Green, minus the connections to a twelfth input. 216 { 217 s.conditionalSwap!(less, Range, 218 0,1, 2,3, 4,5, 6,7, 8,9, 219 0,2, 1,3, 4,6, 5,7, 8,10, 220 1,2, 5,6, 9,10, 221 1,5, 6,10, 222 2,6, 5,9, 223 0,4, 1,5, 6,10, 224 3,7, 4,8, 225 0,4, 226 1,4, 3,8, 7,10, 227 2,3, 8,9, 228 2,4, 3,5, 6,8, 7,9, 229 3,4, 5,6, 7,8); 230 } 231 else static if (n == 12) // Shapiro and Green. 232 { 233 s.conditionalSwap!(less, Range, 234 0,1, 2,3, 4,5, 6,7, 8,9, 10,11, 235 0,2, 1,3, 4,6, 5,7, 9,11, 8,10, 236 1,2, 5,6, 9,10, 237 1,5, 6,10, 238 2,6, 5,9, 239 0,4, 1,5, 6,10, 7,11, 240 3,7, 4,8, 241 0,4, 7,11, 242 1,4, 3,8, 7,10, 243 2,3, 8,9, 244 2,4, 3,5, 6,8, 7,9, 245 3,4, 5,6, 7,8); 246 } 247 else static if (n == 13) // Generated by the END algorithm. 248 { 249 s.conditionalSwap!(less, Range, 250 0,12, 1,7, 2,6, 3,4, 5,8, 9,11, 251 0,1, 2,3, 4,6, 5,9, 7,12, 8,11, 252 0,2, 1,4, 3,7, 6,12, 7,8, 10,11, 253 4,9, 6,10, 11,12, 254 1,7, 3,4, 5,6, 8,9, 10,11, 255 1,3, 2,6, 4,7, 8,10, 9,11, 256 0,5, 2,5, 6,8, 9,10, 257 1,2, 3,5, 4,6, 7,8, 258 2,3, 4,5, 6,7, 8,9, 259 3,4, 5,6); 260 } 261 else static if (n == 14) // Green's construction for 16 inputs minus connections to the fifteenth and sixteenth inputs. 262 { 263 s.conditionalSwap!(less, Range, 264 0,1, 2,3, 4,5, 6,7, 8,9, 10,11, 12,13, 265 0,2, 1,3, 4,6, 5,7, 8,10, 9,11, 266 0,4, 1,5, 2,6, 3,7, 8,12, 9,13, 267 0,8, 1,9, 2,10, 3,11, 4,12, 5,13, 268 3,12, 5,10, 6,9, 269 1,2, 4,8, 7,11, 270 1,4, 2,8, 7,13, 271 2,4, 3,8, 5,6, 7,12, 9,10, 11,13, 272 3,5, 6,8, 7,9, 10,12, 273 3,4, 5,6, 7,8, 9,10, 11,12, 274 6,7, 8,9); 275 } 276 else static if (n == 15) // Green's construction for 16 inputs minus connections to the sixteenth input. 277 { 278 s.conditionalSwap!(less, Range, 279 0,1, 2,3, 4,5, 6,7, 8,9, 10,11, 12,13, 280 0,2, 1,3, 4,6, 5,7, 8,10, 9,11, 12,14, 281 0,4, 1,5, 2,6, 3,7, 8,12, 9,13, 10,14, 282 0,8, 1,9, 2,10, 3,11, 4,12, 5,13, 6,14, 283 3,12, 5,10, 6,9, 13,14, 284 1,2, 4,8, 7,11, 285 1,4, 2,8, 7,13, 11,14, 286 2,4, 3,8, 5,6, 7,12, 9,10, 11,13, 287 3,5, 6,8, 7,9, 10,12, 288 3,4, 5,6, 7,8, 9,10, 11,12, 289 6,7, 8,9); 290 } 291 else static if (n == 16) // Green's construction. TODO Use 10-step Bitonic sorter instead? 292 { 293 s.conditionalSwap!(less, Range, 294 0,1, 2,3, 4,5, 6,7, 8,9, 10,11, 12,13, 14,15, 295 0,2, 1,3, 4,6, 5,7, 8,10, 9,11, 12,14, 13,15, 296 0,4, 1,5, 2,6, 3,7, 8,12, 9,13, 10,14, 11,15, 297 0,8, 1,9, 2,10, 3,11, 4,12, 5,13, 6,14, 7,15, 298 3,12, 5,10, 6,9, 13,14, 299 1,2, 4,8, 7,11, 300 1,4, 2,8, 7,13, 11,14, 301 2,4, 3,8, 5,6, 7,12, 9,10, 11,13, 302 3,5, 6,8, 7,9, 10,12, 303 3,4, 5,6, 7,8, 9,10, 11,12, 304 6,7, 8,9); 305 } 306 else static if (n == 18) // Baddar's PHD thesis, chapter 6. Fewest stages but 2 comparators more than 'batcher' 307 { 308 s.conditionalSwap!(less, Range, 309 0,1, 2,3, 4,5, 6,7, 8,9, 10,11, 12,13, 14,15, 16,17, 310 0,2, 1,3, 4,6, 5,7, 8,10, 9,11, 12,17, 13,14, 15,16, 311 0,4, 1,5, 2,6, 3,7, 9,10, 8,12, 11,16, 13,15, 14,17, 312 0,8, 1,13, 2,4, 3,5, 6,17, 7,16, 9,15, 10,14, 11,12, 313 0,1, 2,11, 3,15, 4,10, 5,12, 6,13, 7,14, 8,9, 16,17, 314 1,8, 3,10, 4,15, 5,11, 6,9, 7,13, 14,16, 315 1,2, 3,6, 4,8, 5,9, 7,11, 10,13, 12,16, 14,15, 316 2,3, 5,8, 6,7, 9,10, 11,14, 12,13, 317 2,4, 3,5, 6,8, 7,9, 10,11, 12,14, 13,15, 318 3,4, 5,6, 7,8, 9,10, 11,12, 13,14, 319 4,5, 6,7, 8,9, 10,11, 12,13); 320 } 321 else static if (n == 22) // Baddar's PHD thesis, chapter 7. Fewest stages but 2 comparators more than 'batcher' 322 { 323 s.conditionalSwap!(less, Range, 324 0,1, 2,3, 4,5, 6,7, 8,9, 10,11, 12,13, 14,15, 16,17, 18,19, 20,21, 325 0,5, 1,3, 2,4, 6,8, 7,9, 10,12, 11,13, 14,16, 15,17, 18,20, 19,21, 326 6,10, 7,11, 8,12, 9,13, 14,18, 15,19, 16,20, 17,21, 327 0,2, 1,4, 3,5, 7,15, 8,16, 9,17, 11,19, 328 0,10, 1,18, 3,12, 5,20, 13,21, 329 0,7, 2,4, 3,15, 6,14, 9,18, 17,20, 330 0,6, 1,8, 2,11, 3,8, 4,16, 5,10, 12,19, 13,14, 20,21, 331 2,13, 4,7, 5,9, 10,15, 11,17, 12,18, 14,16, 16,20, 18,19, 332 1,3, 2,6, 4,5, 7,9, 8,13, 10,11, 12,14, 15,17, 19,20, 333 3,5, 4,6, 7,8, 9,13, 10,12, 11,14, 15,18, 16,17, 334 1,2, 5,10, 6,7, 8,9, 11,12, 13,15, 14,16, 18,19, 335 2,3, 5,7, 8,10, 9,11, 12,13, 14,15, 16,18, 17,19, 336 2,4, 3,6, 7,8, 9,10, 11,12, 13,14, 15,16, 17,18, 337 3,4, 5,6, 8,9, 10,11, 12,13, 14,15, 16,17, 338 4,5, 6,7); 339 } 340 else 341 { 342 static assert(0, "Unsupported length: " ~ n.stringof); 343 } 344 345 import std.range : assumeSorted; 346 return s.assumeSorted!less; 347 } 348 349 /** Sort range `x` of length `n` using a networking sort. 350 */ 351 auto networkSortExactly(uint n, alias less = "a < b", Range)(Range r) 352 if (n >= 2 && 353 isRandomAccessRange!Range) 354 in 355 { 356 assert(r.length == n); 357 } 358 do 359 { 360 x[].networkSortUpTo!(n); 361 } 362 363 /** Sort static array `x` of length `n` using a networking sort. 364 */ 365 auto networkSortExactly(alias less = "a < b", T, size_t n)(ref T[n] x) 366 if (n >= 2) 367 { 368 x[].networkSortUpTo!(n); 369 } 370 371 /// 372 @safe pure nothrow @nogc unittest 373 { 374 int[4] x = [2, 3, 0, 1]; 375 const int[4] y = [0, 1, 2, 3]; 376 x.networkSortExactly; 377 assert(x[].equal(y[])); 378 } 379 380 /** Hybrid sort `r` using `networkSortUpTo` if length of `r` is 381 * less-than-or-equal to `networkSortMaxLength` and `std.algorithm.sorting.sort` 382 * otherwise. 383 */ 384 auto hybridSort(alias less = "a < b", Range)(Range r) 385 if (isRandomAccessRange!Range) 386 { 387 import std.algorithm.sorting : isSorted; 388 static foreach (n; 2 .. networkSortMaxLength + 1) 389 { 390 static if (__traits(compiles, { r.networkSortUpTo!(n, less); })) 391 { 392 if (n == r.length) 393 { 394 auto s = r.networkSortUpTo!(n, less); 395 assert(s.isSorted!less); 396 return s; 397 } 398 } 399 } 400 import std.algorithm.sorting : sort; 401 return sort!less(r); 402 } 403 404 /// 405 @safe pure unittest 406 { 407 import std.algorithm.sorting : isSorted; 408 import std.algorithm.iteration : permutations; 409 import std.range : iota; 410 import std.random : randomShuffle, Random; 411 412 Random random; 413 414 alias T = uint; 415 416 const maxFullPermutationTestLength = 8; 417 const maxTriedShufflings = 10_000; // maximum number of shufflings to try 418 419 import std.meta : AliasSeq; 420 foreach (less; AliasSeq!("a < b", "a > b")) 421 { 422 foreach (const n; 0 .. networkSortMaxLength + 1) 423 { 424 if (n > maxFullPermutationTestLength) // if number of elements is too large 425 { 426 foreach (x; 0 .. maxTriedShufflings) 427 { 428 import std.array : array; 429 auto y = iota(0, n).array; 430 y.randomShuffle(random); 431 y.hybridSort!less; 432 assert(y.isSorted!less); 433 } 434 } 435 else 436 { 437 foreach (x; iota(0, n).permutations) 438 { 439 import std.array : array; 440 auto y = x.array; 441 y.hybridSort!less; 442 assert(y.isSorted!less); 443 } 444 } 445 } 446 } 447 }