1(**** ML Programs from Chapter 3 of
2
3  ML for the Working Programmer, 2nd edition
4  by Lawrence C. Paulson, Computer Laboratory, University of Cambridge.
5  (Cambridge University Press, 1996)
6
7Copyright (C) 1996 by Cambridge University Press.
8Permission to copy without fee is granted provided that this copyright
9notice and the DISCLAIMER OF WARRANTY are included in any copy.
10
11DISCLAIMER OF WARRANTY.  These programs are provided `as is' without
12warranty of any kind.  We make no warranties, express or implied, that the
13programs are free of error, or are consistent with any particular standard
14of merchantability, or that they will meet your requirements for any
15particular application.  They should not be relied upon for solving a
16problem whose incorrect solution could result in injury to a person or loss
17of property.  If you do use the programs or functions in such a manner, it
18is at your own risk.  The author and publisher disclaim all liability for
19direct, incidental or consequential damages resulting from your use of
20these programs or functions.
21****)
22
23
24(*** Making change.  Thanks to Dr Martin Richard for this idea! ***)
25
26(*Return first way of making change; NOT EXHAUSTIVE*)
27fun change (coinvals, 0)         = []
28  | change (c::coinvals, amount) =
29      if amount<c then change(coinvals, amount)
30                  else c :: change(c::coinvals, amount-c);
31
32val gb_coins = [50,20,10,5,2,1]
33and us_coins = [25,10,5,1];
34
35(*Return all ways of making change (not just the first way)  *)
36fun allChange (coins, coinvals, 0)         = [coins]
37  | allChange (coins, [],    amount)       = []
38  | allChange (coins, c::coinvals, amount) =
39      if amount<0 then []
40      else allChange(c::coins, c::coinvals, amount-c) @
41	   allChange(coins, coinvals, amount);
42
43
44(*** Binary Arithmetic ***)
45
46fun bincarry (0, ps) = ps
47  | bincarry (1, []) = [1]
48  | bincarry (1, p::ps) = (1-p) :: bincarry(p, ps);
49
50(*sum of two binary numbers with carry*)
51fun binsum (c, [], qs) = bincarry (c,qs)
52  | binsum (c, ps, []) = bincarry (c,ps)
53  | binsum (c, p::ps, q::qs) =
54      ((c+p+q) mod 2) :: binsum((c+p+q) div 2, ps, qs);
55
56(*product of two binary numbers*)
57fun binprod ([], _) = []
58  | binprod (0::ps, qs) = 0::binprod(ps,qs)
59  | binprod (1::ps, qs) = binsum(0, qs, 0::binprod(ps,qs));
60
61
62(*** Matrix transpose ***)
63
64fun headcol []    = []
65  | headcol ((x::_) :: rows) = x :: headcol rows;
66
67fun tailcols []    = []
68  | tailcols ((_::xs) :: rows) = xs :: tailcols rows;
69
70(*LOOPS given empty list!*)
71fun transp ([]::rows) = []
72  | transp rows = headcol rows :: transp (tailcols rows);
73
74
75(*** Matrix product ***)
76
77fun dotprod([], []) = 0.0
78  | dotprod(x::xs,y::ys) = x*y + dotprod(xs,ys);
79
80fun rowprod(row, []) = []
81  | rowprod(row, col::cols) =
82        dotprod(row,col) :: rowprod(row,cols);
83
84fun rowlistprod([], cols) = []
85  | rowlistprod(row::rows, cols) =
86        rowprod(row,cols) :: rowlistprod(rows,cols);
87
88fun matprod(rowsA,rowsB) = rowlistprod(rowsA, transp rowsB);
89
90
91(*** Gaussian Elimination (from Sedgewick, Chapter 5) ***)
92
93(*the first row with absolutely greatest head*)
94fun pivotrow [row] = row : real list
95  | pivotrow (row1::row2::rows) =
96      if abs(hd row1) >= abs(hd row2)
97      then pivotrow(row1::rows)
98      else pivotrow(row2::rows);
99
100(*the matrix excluding the first row with head p*)
101fun delrow (p, []) = []
102  | delrow (p, row::rows) =
103      if Real.==(p, hd row) then rows
104      else row :: delrow(p, rows);
105
106fun scalarprod(k, []) = [] : real list
107  | scalarprod(k, x::xs) = k*x :: scalarprod(k,xs);
108
109fun vectorsum ([], []) = [] : real list
110  | vectorsum (x::xs,y::ys) = x+y :: vectorsum(xs,ys);
111
112fun gausselim [row] = [row]
113  | gausselim rows =
114      let val p::prow = pivotrow rows
115          fun elimcol [] = []
116            | elimcol ((x::xs)::rows) =
117                  vectorsum(xs, scalarprod(~x/p, prow))
118                  :: elimcol rows
119      in  (p::prow) :: gausselim(elimcol(delrow(p,rows)))
120      end;
121
122fun solutions [] = [~1.0]
123  | solutions((x::xs)::rows) =
124      let val solns = solutions rows
125      in ~(dotprod(solns,xs)/x) :: solns  end;
126
127
128(*** Dijkstra's problems ***)
129
130(** Writing a number as the sum of two squares.
131    A number is the sum of two squares iff
132    its square-free odd prime factors all are congruent to 1 (mod 4) --
133    thus 2, 5, 13, 17, 29, 37, 41,  53,  61, 73, 89, 97, 101, 109, 113...
134    See H. Davenport (1952), chapter V.**)
135
136fun squares r =
137  let fun between (x,y) =  (*all pairs between x and y*)
138        let val diff = r - x*x
139            fun above y =  (*all pairs above y*)
140                if y>x then []
141                else if y*y<diff then above (y+1)
142                else if y*y=diff then (x,y)::between(x-1,y+1)
143                else (* y*y>diff *)  between(x-1,y)
144        in above y  end;
145      val firstx = floor(Math.sqrt(real r))
146  in between (firstx, 0) end;
147
148(** the next permutation **)
149
150fun next(xlist, y::ys) : int list =
151    if hd xlist <= y then  next(y::xlist, ys)  (*still increasing*)
152    else  (*swap y with greatest xk such that x>=xk>y *)
153      let fun swap [x] = y::x::ys
154            | swap (x::xk::xs) =          (*x >= xk *)
155                if xk>y then x::swap(xk::xs)
156                else (y::xk::xs)@(x::ys)
157                         (* x > y >= xk >= xs *)
158      in swap(xlist) end;
159
160fun nextperm (y::ys) = next([y], ys);
161
162
163(*** Finite sets ***)
164
165(*membership in a list*)
166infix mem;
167fun x mem []  =  false
168  | x mem (y::l)  =  (x=y) orelse (x mem l);
169
170(*insertion into list if not already there*)
171fun newmem(x,xs) = if x mem xs then  xs   else  x::xs;
172
173(*insert the list of xs into the ys, adding no duplicates*)
174fun union([],ys) = ys
175  | union(x::xs, ys) = newmem(x, union(xs, ys));
176
177fun inter([],ys) = []
178  | inter(x::xs, ys) =
179        if x mem ys then x::inter(xs, ys)
180                    else    inter(xs, ys);
181
182fun powset ([], base) = [base]
183  | powset (x::xs, base) = powset(xs, base) @ powset(xs, x::base);
184
185fun cartprod ([], ys) = []
186  | cartprod (x::xs, ys) =
187        let val rest = cartprod(xs,ys)
188            fun pairx [] = rest
189              | pairx(y::ytail) = (x,y) :: (pairx ytail)
190        in  pairx ys  end;
191
192
193(*** Graph algorithms ***)
194
195fun nexts (a, []) = []
196  | nexts (a, (x,y)::pairs) =
197      if a=x then  y :: nexts(a,pairs)
198             else       nexts(a,pairs);
199
200fun depthf ([], graph, visited) = rev visited
201  | depthf (x::xs, graph, visited) =
202      if x mem visited then depthf (xs, graph, visited)
203      else depthf (nexts(x,graph) @ xs, graph, x::visited);
204
205(*Alternative, faster function for depth-first search*)
206fun depth ([], graph, visited) = rev visited
207  | depth (x::xs, graph, visited) =
208      depth (xs, graph,
209             if x mem visited  then  visited
210             else depth (nexts(x,graph), graph, x::visited));
211
212fun topsort graph =
213  let fun sort ([], visited) = visited
214        | sort (x::xs, visited) =
215            sort(xs, if x mem visited  then  visited
216                     else x :: sort(nexts(x,graph), visited));
217      val (xs,_) = ListPair.unzip graph
218  in sort(xs, []) end;
219
220fun pathsort graph =
221  let fun sort ([], path, visited) = visited
222        | sort (x::xs, path, visited) =
223            if x mem path then hd[] (*abort!!*)
224            else sort(xs, path,
225                      if x mem visited  then  visited
226                      else x :: sort(nexts(x,graph), x::path, visited))
227      val (xs,_) = ListPair.unzip graph
228  in sort(xs, [], []) end;
229
230
231fun newvisit (x, (visited,cys)) = (x::visited, cys);
232
233fun cyclesort graph =
234  let fun sort ([], path, (visited,cys)) = (visited, cys)
235        | sort (x::xs, path, (visited,cys)) =
236            sort(xs, path,
237               if x mem path   then  (visited, x::cys)
238               else if x mem visited then (visited, cys)
239               else newvisit(x, sort(nexts(x,graph),
240                                     x::path, (visited,cys))))
241      val (xs,_) = ListPair.unzip graph
242  in sort(xs, [], ([],[])) end;
243
244
245(*** Sorting ***)
246
247(** Random numbers, courtesy Stephen K. Park and Keith W. Miller,
248	CACM 31 (1988), 1192-1201.  **)
249local val a = 16807.0  and  m = 2147483647.0
250in  fun nextrandom seed =
251          let val t = a*seed
252          in  t - m * real(floor(t/m))  end
253
254    (*truncate to integer from 1 to k*)
255    and truncto k r = 1 + floor((r / m) * (real k))
256end;
257
258fun randlist (n,seed,seeds) =
259    if n=0  then  (seed,seeds)
260    else  randlist(n-1, nextrandom seed, seed::seeds);
261
262
263(** insertion sort: non-iterative is faster **)
264fun ins (x, []): real list = [x]
265  | ins (x, y::ys) =
266      if x<=y then x::y::ys    (*it belongs here*)
267              else y::ins(x,ys);
268
269fun insort [] = []
270  | insort (x::xs) = ins(x, insort xs);
271
272
273(*quicksort*)
274fun quick [] = []
275  | quick [x] = [x]
276  | quick (a::bs) =  (*the head "a" is the pivot*)
277      let fun partition (left,right,[]) : real list =
278                (quick left) @ (a :: quick right)
279            | partition (left,right, x::xs) =
280                if x<=a then partition (x::left, right, xs)
281                        else partition (left, x::right, xs)
282      in  partition([],[],bs)  end;
283
284
285(** Top-down merge sort **)
286
287fun merge([],ys) = ys : real list
288  | merge(xs,[]) = xs
289  | merge(x::xs, y::ys) =
290      if x<=y then x::merge(xs,  y::ys)
291              else y::merge(x::xs,  ys);
292
293(*naive version -- like Bird and Wadler, following Sedgewick*)
294fun tmergesort [] = []
295  | tmergesort [x] = [x]
296  | tmergesort xs =
297      let val k = length xs div 2
298      in  merge (tmergesort (List.take(xs,k)),
299                 tmergesort (List.drop(xs,k)))
300      end;
301
302(*faster version*)
303fun tmergesort' xs =
304      let fun sort (0, xs) = ([], xs)
305	    | sort (1, x::xs) = ([x], xs)
306	    | sort (n, xs) =
307		let val (l1, xs1) = sort ((n+1) div 2, xs)
308		    val (l2, xs2) = sort (n div 2, xs1)
309		in (merge (l1,l2), xs2)
310		end
311          val (l, _) = sort (length xs, xs)
312      in l end;
313
314
315(** Bottom-up merge sort **)
316
317fun mergepairs([l], k) = [l]
318  | mergepairs(l1::l2::ls, k) =
319      if k mod 2 = 1 then l1::l2::ls
320      else mergepairs(merge(l1,l2)::ls, k div 2);
321
322fun sorting([], ls, r) = hd(mergepairs(ls,0))
323  | sorting(x::xs, ls, r) = sorting(xs, mergepairs([x]::ls, r+1), r+1);
324
325fun sort xs = sorting(xs, [[]], 0);
326
327(*O'Keefe's samsort*)
328fun nextrun(run, []) =       (rev run, []: real list)
329  | nextrun(run, x::xs) =
330        if  x < hd run then  (rev run, x::xs)
331                       else  nextrun(x::run, xs);
332
333fun samsorting([], ls, k) = hd(mergepairs(ls,0))
334  | samsorting(x::xs, ls, k) =
335      let val (run, tail) = nextrun([x], xs)
336      in  samsorting(tail, mergepairs(run::ls, k+1), k+1)
337      end;
338
339fun samsort xs = samsorting(xs, [[]], 0);
340
341
342(*** Polynomial arithmetic -- as a structure ***)
343
344(*An obscure built-in infix in some compilers; we need it nonfix*)
345nonfix rem;
346
347structure Poly =
348  struct
349  type t = (int*real)list;
350
351  val zero = [];
352
353  (** Sum of two polynomials **)
354  fun sum ([], us)               = us : t
355    | sum (ts, [])               = ts
356    | sum ((m,a)::ts, (n,b)::us) =
357	     if m>n then (m,a) :: sum (ts, (n,b)::us)
358	else if n>m then (n,b) :: sum (us, (m,a)::ts)
359	else (*m=n*)
360	     if Real.==(a+b, 0.0) then sum (ts, us)
361			          else (m, a+b) :: sum (ts, us);
362
363  (** Product of a term and a polynomial **)
364  fun termprod ((m,a), [])        = [] : t
365    | termprod ((m,a), (n,b)::ts) =
366	(m+n, a*b) :: termprod ((m,a), ts);
367
368  (** Product of two polynomials **)
369
370  (*Product by balanced merging; could improve speed, like with merge sort*)
371  fun prod ([], us)      = []
372    | prod ([(m,a)], us) = termprod ((m,a), us)
373    | prod (ts, us)      =
374	let val k = length ts div 2
375	in  sum (prod (List.take(ts,k), us), prod (List.drop(ts,k), us))
376	end;
377
378  (** Division: quotient and remainder **)
379
380  (*Division by zero -- empty polynomial -- raises exception Match*)
381  fun quorem (ts, (n,b)::us) =
382    let fun dividing ([],        qs) = (rev qs, [])
383	  | dividing ((m,a)::ts, qs) =
384	      if m<n then (rev qs, (m,a)::ts)
385	      else dividing (sum (ts, termprod ((m-n, ~a/b), us)),
386			(m-n, a/b) :: qs)
387    in  dividing (ts, [])  end;
388
389  fun quo (ts,us) = #1(quorem (ts,us))
390  and rem (ts,us) = #2(quorem (ts,us));
391
392  (*A bad gcd algorithm: results are not in primitive form*)
393  fun gcd ([], us) = us
394    | gcd (ts, us) = gcd (rem (us,ts), ts);
395
396  end;
397
398
399