sorting.lisp 11.7 KB
Newer Older
1 2
;; Sorting
;; Liam Healy, Fri Apr 14 2006 - 20:20
3
;; Time-stamp: <2012-01-13 12:01:10EST sorting.lisp>
4
;;
5
;; Copyright 2006, 2007, 2008, 2009, 2011, 2012 Liam M. Healy
6 7 8 9 10 11 12 13 14 15 16 17 18 19
;; Distributed under the terms of the GNU General Public License
;;
;; This program is free software: you can redistribute it and/or modify
;; it under the terms of the GNU General Public License as published by
;; the Free Software Foundation, either version 3 of the License, or
;; (at your option) any later version.
;;
;; This program is distributed in the hope that it will be useful,
;; but WITHOUT ANY WARRANTY; without even the implied warranty of
;; MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
;; GNU General Public License for more details.
;;
;; You should have received a copy of the GNU General Public License
;; along with this program.  If not, see <http://www.gnu.org/licenses/>.
liam's avatar
liam committed
20 21 22 23

(in-package :gsl)

;;; #'heapsort has just a cursory port, use CL's #'sort.
24
;;; [2012-01-03 Tue 13:15] To do: rewrite other index functions to take an optional array or size (sort-smallest-index, sort-largest-index already done).
liam's avatar
liam committed
25 26 27 28 29 30 31 32 33

;;;;****************************************************************************
;;;; Heapsort, not recommended
;;;;****************************************************************************

(defmacro defcomparison (name &body body)
  `(cffi:defcallback ,name :int ((a :pointer) (b :pointer))
    ,body))

liam's avatar
liam committed
34
(defmfun heapsort (array count size function)
35
  "gsl_heapsort"
36
  ((array :pointer) (count :sizet) (size :sizet) (function :pointer))
37
  :documentation			; FDL
38 39 40
  "Sort the count elements of the array of size specified
   into ascending order using the comparison
   function.  The type of the comparison function is defined by,
liam's avatar
liam committed
41
   A comparison function should return a negative integer if the first
42
   argument is less than the second argument, zero if the two arguments
liam's avatar
liam committed
43 44
   are equal and a positive integer if the first argument is greater than
   the second argument."
45
  :c-return :void)
liam's avatar
liam committed
46

liam's avatar
liam committed
47
(defmfun heapsort-index (p array count size function)
liam's avatar
liam committed
48
  "gsl_heapsort_index"
49
  ((p :sizet) (array :pointer) (count :sizet) (size :sizet) (function :pointer))
50 51 52 53 54
  :documentation			; FDL
  "Indirectly sort the count elements of the array
   array, each of size given, into ascending order using the
   comparison function.  The resulting permutation is stored
   in p, an array of length n.  The elements of p give the
liam's avatar
liam committed
55
   index of the array element which would have been stored in that position
56 57 58 59
   if the array had been sorted in place.  The first element of p
   gives the index of the least element in array, and the last
   element of p gives the index of the greatest element in
   array.  The array itself is not changed.")
liam's avatar
liam committed
60 61

;;;;****************************************************************************
62
;;;; Sorting vectors
liam's avatar
liam committed
63 64
;;;;****************************************************************************

65 66
;;; The GSL vector sorting routines come in two forms, those that sort
;;; on GSL vectors (with "vector" in the name), and those that sort on
67 68 69 70 71
;;; C arrays.  The former can sort vectors, the latter can both
;;; vectors and matrices.  The functions #'sort-index,
;;; #'sort-smallest-index, #'sort-largest-index could be made to work
;;; on matrices, but the indexing would have to be worked out
;;; correctly.
liam's avatar
liam committed
72

73 74 75
;;; It ought to be possible to provide a stride argument, but this
;;; gives an error:
;;;(defmfun msort ((v vector) &optional (stride 1))
76

77
(defmfun msort ((v both))
78
  ("gsl_sort" :type)
79
  (((grid:foreign-pointer v) :pointer) (1 :sizet) ((size v) :sizet))
80 81 82
  :definition :generic
  :element-types :no-complex
  :c-return :void
83
  :inputs (v)
84 85 86 87 88 89 90 91 92 93 94
  :outputs (v)
  :documentation			; FDL
  "Sort the n elements of the array data with stride stride into
   ascending numerical order.")

(defmfun sort-vector ((v vector))
  ("gsl_sort_vector" :type)
  (((mpointer v) :pointer))
  :definition :generic
  :element-types :no-complex
  :c-return :void
95
  :inputs (v)
96 97 98 99
  :outputs (v)
  :documentation			; FDL
  "Sort the elements of the vector v into ascending numerical order.")

100
(defmfun sort-index ((permutation permutation) (vector vector))
101 102
  ("gsl_sort" :type "_index")
  (((mpointer permutation) :pointer)
103
   ((grid:foreign-pointer vector) :pointer)
104
   (1 :sizet) ((dim0 vector) :sizet))
105 106 107 108 109 110 111 112 113 114 115 116
  :definition :generic
  :element-types :no-complex
  :c-return :void
  :inputs (vector)
  :outputs (permutation)
  :documentation			; FDL
  "Indirectly sort the n elements of the array vector with stride stride
   into ascending order, storing the resulting permutation.  The latter
   must be created with the same size as the vector.
   The elements of permutation give the index of the
   array element which would have been stored in that position if the
   array had been sorted in place. The array data is not changed.")
liam's avatar
liam committed
117

118 119 120 121 122 123 124 125 126 127 128
(defmfun sort-vector-index ((permutation permutation) (vector vector))
  ("gsl_sort_vector" :type "_index")
  (((mpointer permutation) :pointer)
   ((mpointer vector) :pointer))
  :definition :generic
  :element-types :no-complex
  :c-return :void
  :inputs (vector)
  :outputs (permutation)
  :documentation			; FDL
  "Indirectly sort the elements of the vector v into
129 130 131 132 133
   ascending order, storing the resulting permutation in p.  The
   elements of p give the index of the vector element which would
   have been stored in that position if the vector had been sorted in
   place.  The first element of p gives the index of the least element
   in v and the last element of p gives the index of the
134 135 136 137 138 139 140 141
   greatest element in v.  The vector v is not changed.")

;;;;****************************************************************************
;;;; Selecting the k smallest or largest elements
;;;;****************************************************************************

(defmfun sort-vector-smallest (dest (v vector))
  ("gsl_sort_vector" :type "_smallest")
142
  (((grid:foreign-pointer dest) :pointer) ((dim0 dest) :sizet)
143 144 145 146
   ((mpointer v) :pointer))
  :definition :generic
  :element-types :no-complex
  :c-return :void
147
  :inputs (v)
148 149 150 151 152
  :outputs (dest)
  :documentation			; FDL
  "Find the smallest elements of the vector v and put them into dest,
   which must be shorter than v.")

153
(defmfun sort-smallest (dest (v both))
154
  ("gsl_sort" :type "_smallest")
155
  (((grid:foreign-pointer dest) :pointer) ((size dest) :sizet)
156
   ((grid:foreign-pointer v) :pointer)
157 158
   (1 :sizet)				; stride, set to 1 for now
   ((size v) :sizet))
159 160
  :definition :generic
  :element-types :no-complex
161
  :c-return :void
162
  :inputs (v)
163 164 165 166
  :outputs (dest)
  :documentation			; FDL
  "Find the smallest elements of the vector v and put them into dest,
   which must be shorter than v.")
167

168
(defmfun sort-vector-smallest-index (combination (v vector))
169
  ("gsl_sort_vector" :type "_smallest_index")
170
  (((grid:foreign-pointer combination) :pointer) ((size combination) :sizet)
171
   ((mpointer v) :pointer)
172 173
   (1 :sizet)				; stride, set to 1 for now
   ((first (grid:dimensions combination)) :sizet))
174 175 176
  :definition :generic
  :element-types :no-complex
  :c-return :void
177 178
  :inputs (v)
  :outputs (combination)
179 180
  :documentation
  "The indices of the smallest elements of the vector stored,
181 182 183
    returned as a CL vector of element type fixnum.  If
    indices is a positive initeger, a vector will be
    allocated and returned.  If it is a CL vector,
184 185
    it will be filled with the indices.")

186 187
(defmfun sort-smallest-index
    ((v vector) &optional (size-or-array 3)
188
		&aux (combination (vdf size-or-array ':sizet)))
189
  ("gsl_sort" :type "_smallest_index")
190
  (((grid:foreign-pointer combination) :pointer) ((size combination) :sizet)
191
   ((grid:foreign-pointer v) :pointer)
192 193
   (1 :sizet)				; stride, set to 1 for now
   ((first (grid:dimensions v)) :sizet))
194
  :element-types :no-complex
195
  :definition :generic
196
  :c-return :void
197 198
  :inputs (v)
  :outputs (combination)
199
  :documentation
200
  "The indices of the smallest elements of the vector.  If size-or-array is an integer, it is the number of smallest elements.  If it is an array of :sizet elements, it is filled.")
201 202 203

(defmfun sort-vector-largest (dest (v vector))
  ("gsl_sort_vector" :type "_largest")
204
  (((grid:foreign-pointer dest) :pointer) ((dim0 dest) :sizet)
205 206 207 208
   ((mpointer v) :pointer))
  :definition :generic
  :element-types :no-complex
  :c-return :void
209
  :inputs (v)
210 211 212 213 214
  :outputs (dest)
  :documentation			; FDL
  "Find the largest elements of the vector v and put them into dest,
   which must be shorter than v.")

215
(defmfun sort-largest (dest (v both))
216
  ("gsl_sort" :type "_largest")
217
  (((grid:foreign-pointer dest) :pointer) ((size dest) :sizet)
218
   ((grid:foreign-pointer v) :pointer)
219 220
   (1 :sizet)				; stride, set to 1 for now
   ((size v) :sizet))
221 222
  :definition :generic
  :element-types :no-complex
223
  :c-return :void
224
  :inputs (v)
225 226 227 228 229
  :outputs (dest)
  :documentation			; FDL
  "Find the largest elements of the vector v and put them into dest,
   which must be shorter than v.")

230
(defmfun sort-vector-largest-index (combination (v vector))
231
  ("gsl_sort_vector" :type "_largest_index")
232
  (((grid:foreign-pointer combination) :pointer) ((size combination) :sizet)
233
   ((mpointer v) :pointer)
234 235
   (1 :sizet)				; stride, set to 1 for now
   ((first (grid:dimensions combination)) :sizet))
236 237
  :definition :generic
  :element-types :no-complex
238
  :c-return :void
239 240
  :inputs (v)
  :outputs (combination)
241 242
  :documentation
  "The indices of the largest elements of the vector stored,
243 244 245
    returned as a CL vector of element type fixnum.  If
    indices is a positive initeger, a vector will be
    allocated and returned.  If it is a CL vector,
246
    it will be filled with the indices.")
247

248 249
(defmfun sort-largest-index
    ((v vector) &optional (size-or-array 3)
250
		&aux (combination (vdf size-or-array ':sizet)))
251
  ("gsl_sort" :type "_largest_index")
252
  (((grid:foreign-pointer combination) :pointer) ((size combination) :sizet)
253
   ((grid:foreign-pointer v) :pointer)
254 255
   (1 :sizet)				; stride, set to 1 for now
   ((first (grid:dimensions v)) :sizet))
256 257
  :element-types :no-complex
  :definition :generic
258
  :c-return :void
259 260
  :inputs (v)
  :outputs (combination)
261
  :documentation
262
  "The indices of the largest elements of the vector.  If size-or-array is an integer, it is the number of smallest elements.  If it is an array of sizet elements, it is filled.")
liam's avatar
liam committed
263 264

;;;;****************************************************************************
265
;;;; Examples and unit test
liam's avatar
liam committed
266 267
;;;;****************************************************************************

268
(generate-all-array-tests sort-vector :no-complex
269
 (let ((v1 (array-default 8)))
270
   ;; or you can use msort
271
   (grid:copy-to (sort-vector v1))))
272 273

(generate-all-array-tests sort-matrix :no-complex
274
 (let ((m1 (array-default '(3 3))))
275
   (grid:copy-to (msort m1))))
276 277

(generate-all-array-tests sort-vector-index :no-complex
278
 (let ((perm (make-permutation 8))
279 280
	(v1 (array-default 8)))
   (sort-vector-index perm v1)
281
   (grid:copy-to perm)))
282 283

(generate-all-array-tests sort-vector-smallest :no-complex
284
 (let ((v1 (array-default 8))
285
	(v2 (array-default 3)))
286
   (grid:copy-to (sort-vector-smallest v2 v1))))
287 288

(generate-all-array-tests sort-matrix-smallest :no-complex
289
 (let ((m1 (array-default '(3 3)))
290
	(m2 (array-default '(2 3) t)))
291
   (grid:copy-to (sort-smallest m2 m1))))
292 293

(generate-all-array-tests sort-vector-smallest-index :no-complex
Liam Healy's avatar
Liam Healy committed
294
 (let ((comb (make-combination 8 3 nil))
295
	(v1 (array-default 8)))
296
   (grid:copy-to (sort-vector-smallest-index comb v1))))
297 298

(generate-all-array-tests sort-vector-largest :no-complex
299
 (let ((v1 (array-default 8))
300
	(v2 (array-default 3)))
301
   (grid:copy-to (sort-vector-largest v2 v1))))
302 303

(generate-all-array-tests sort-matrix-largest :no-complex
304
 (let ((m1 (array-default '(3 3)))
305
	(m2 (array-default '(2 3) t)))
306
   (grid:copy-to (sort-largest m2 m1))))
liam's avatar
liam committed
307

308
(generate-all-array-tests sort-vector-largest-index :no-complex
Liam Healy's avatar
Liam Healy committed
309
 (let ((comb (make-combination 8 3 nil))
310
	(v1 (array-default 8)))
311
   (grid:copy-to (sort-vector-largest-index comb v1))))