The Salty Economist

Things I Should Have Learned in High School
posts - 56, comments - 0, trackbacks - 0

F# Heap Sort

Today's addition to my math library is a heap sort function.

I often run into a problem with collections, where I have a collection of class objects that I need to access based on one of their characteristics.  For example, suppose I have a collection of 10 automobile objects and I wan them sorted by their MPG.

My strategy here is this:

(1)  Create an array of floating point numbers (the MPG for each car)

(2)  Create a parallel array of integers that represent each car's index in the collection.

(3)  Sort the floating point array by MPG, but use the auxiliary integer array to store the ordering.  This allows me to have access the collection (by index) according the the sorted integer array.

So, the first thing we do is set up our car list:

type auto = {model : string; year : int; color : string; mpg : float }

let cars = [{model = "Mustang"; year = 1969; color = "Black"; mpg = 11.0}; 
            {model = "GTO"; year = 1966; color = "Primer Gray"; mpg = 8.5};
            {model = "Cougar"; year = 1968; color = "Red"; mpg = 9.5};
            {model = "VW Bug"; year = 1971; color = "Lime"; mpg = 28.0};
            {model = "Camry"; year = 2004; color = "Tan"; mpg = 25.5};
            {model = "Nova"; year = 1974; color = "White"; mpg = 15.0};
            {model = "Charger"; year = 1971; color = "Red"; mpg = 12.0};
            {model = "Dart"; year = 1963; color = "Light Green"; mpg = 11.0};
            {model = "Comet"; year = 1964; color = "Black"; mpg = 9.0};
            {model = "Mark IV"; year = 1970; color = "Silver"; mpg = 7.5} ]

Then we set up our two arrays:

let intArray = [|1..cars.Length|]

let zArray = Array.create cars.Length 0.0

let rec fillArray aList (rArray : float array) (n : int) =
    match aList with
    | h::t -> rArray.[n] <- h.mpg
              fillArray t rArray (n+1)
    | [] -> rArray

let realArray = fillArray cars zArray 0

I first create an array of integers, that goes from 1 to the length of my cars array.  I then create a float array of zeros of the same length.  I then write a little function that iterates over the car list and fills up the real array with the mpg of each car.  I could fill the array using a simple loop, but I wanted to practice using the pattern matching syntax with the 'cons' pattern.

The statement:

    | h::t -> rArray.[n] <- h.mpg
              fillArray t rArray (n+1)

splits the list between its first item (head or h) and the rest of the list (tail or t).  I then fill the nth spot in the array with the head item's mpg.  I then recursively call the fillArray function, passing it the tail and incrementing the array index (n) by one.

We are now ready to write the Heap Sort.

The way a Heap Sort works is that we first build a Heap.  So what is a Heap?  Aside from a pile of junk in the backyard, I mean.  Technically a heap is a hierarchical tree structure in which each parent node has a value greater than either of it's child node(s).

The way a Heap works is that you organize an array that are organized in a way that resembles a binary tree.  The array element with highest value sits at the first position (position 0) of the array.  Positions 1 and 2 are then two elements that are lower in value than the first element.  At positions 3 and 4 are array elements that are lower in value than element 1;  positions 5 and 6 hold two elements that are lower in value than element 2.

The following picture shows a binary heap, where I have noted the array position followed by the array value.  This is the heap structure we want to build from a car data.  That is. position 0 has a value of 28.0, position 1 has a value of 25.5, and so on.

 

Note that in the tree there is no guarantee that element 4 is lower than element 2, even though it is lower in the table.  For example, element 4 could be 20.0 (still less than its parent 25.5) which is greater than element 2 (15.5), but the required relationships of the tree are still intact.

Here is the fundamental key to the Tree:  If we take any element (call it the ith element), we know that it's parent is as position: (i - 1)/2 (round down).  For example, if i = 8, then its parent is 3; if i = 6, its parent is 2.  Likewise, we know that its children are at positions: (i  * 2) + 1 and (i * 2) + 2

We use these rules to reorganize our array.  Suppose we take any value, we can compare it to the values of where its children reside.  If we find that the parent is smaller than either child, we swap the parent for the largest child.  This swap, means that for that node the required parent/child relationship will hold true.

Here is the code to build the Heap:

---------------------------------------------------

let printA (i : int) (rA : float array) =
    printfn "The Heap: %i : %10.1f %10.1f %10.1f %10.1f %10.1f %10.1f %10.1f %10.1f %10.1f %10.1f "
            i rA.[0] rA.[1] rA.[2] rA.[3] rA.[4] rA.[5] rA.[6] rA.[7] rA.[8] rA.[9]

//
// Left Node
//
let left (i : int) =
    i*2 + 1

//
// Right Node
//
let right (i : int) =
    i*2 + 2

//
//Test if Left Node > ith Node
//
let iBig (rA : float array) (iLeft : int) (i : int) (n : int) =
    if iLeft < n then
        if rA.[iLeft] > rA.[i] then
            iLeft
        else
            i
    else
        i

//
//Swap Values
//
let swap (rA : float array) (i : int) (j : int) =
    let x = rA.[i]
    rA.[i] <- rA.[j]
    rA.[j] <- x
    rA

//
//  Recursive function to build the heap
//
let rec BuildTheHeap (rA : float array) (n : int) (i : int) =
    // Recursive loop that works down to zero (including zero)
    match i with
    | -1 -> rA
    | _ -> let rA = heapIt rA i n
           printA i rA
           BuildTheHeap rA  n (i-1)

---------------------------------------------------

So, how does this work?

 The last function 'BuildTheHeap' takes as arguments the array rA, the array length n, and a position in the array i.  In our example, the array length is 10.  The function is initially called at the midpoint of the array (the value 4 in our case).  The function is recursively called, where the position i is decremented until it hit -1.  We do want to include i of zero (0) in our loop.  In effect, this function starts at the midpoint and moves one position to the left at each iteration.

The function 'heapIt' takes a array position i as input and compares it against the values of its two children two children.  What position's do the children hold?   Form above, we know that its children are at positions: (i  * 2) + 1 and (i * 2) + 2.  Hence, we have our two functions 'left' and 'right' that simply return the left and right child nodes of the ith node.  From this, you can see why we start at the midpoint:  any nodes to the right have the midpoint have NO children.

The function iBig compares the ith node to the Left node and returns the relevant index.  Heapit has another function iBigger that tests the right node against iBig and returns the index of the larger value.

Then, we check to see if iBigger not equal to i.  If iBig equals i, it means that the value at i exceeds each of the values of the children, meaning the required relationships for the binary tree hold -- there is nothing to do.  If iBigger does not equal i, it means the value at iBigger is greater than at i.  We need to swap the values in order to reestablish the required binary tree relationship.  This swapping moves the greater value to the left and the lesser value to the right. 

We are not quite done yet.  For example, let's say the value at position 5 is greater than the value at position 2.  We swap the values, thereby demoting 2.  Now, we want to further check to see if the new value at 5 exceeds the values of its children.  If it does, we stop.  But if not, we swap again and reevaluate.  This process continues to the right until we run out of children or the value finds two children with lesser values.

Note:  the printA function just prints the Array to the console.

In agonizing detail, here is how the Heap is built.

This is the initial array:

Position:        0     1     2     3     4     5     6     7     8     9    10
Start:        11.5   8.5   9.5  28.0  25.5  15.0  12.0  11.0   9.0   7.5    <>

The algorithm works on this array (in place) to build the Heap.

Step 1:  Find the midpoint.  Divide the number of elements (less 1) by 2.  We need to subtract 1 because the array is zero based. so that 9 is really the top position not 10.  We then round down. I this value i: 

Step 2:  Starting at the midpoint (position 4), compare the value with the values at a pair of values that are at positions:

left = i + (2*i) + 1
right = i + (2*i) + 2
(positions 9 and 10 in our case)

If the value at position 4 is lower that the value at positions 9 or 10, swap it with the one that is greater.   (In our example position 10 is outside the array boundary, so it is skipped.)   In the example above, 25.5 is greater than 7.5, so we do nothing.

Start:        11.5   8.5   9.5  28.0  25.5  15.0  12.0  11.0   9.0   7.5

We move to the left one position and test the values 28.0 vs 11.0 and 9.0.  Again we have nothing to do.

Start:        11.5   8.5   9.5  28.0  25.5  15.0  12.0  11.0   9.0   7.5

We again move to the left one position and test the values.  Here we have 9.5 versus the pair 15.0 and 12.0.

The Heap:     11.5   8.5  15.0  28.0  25.5   9.5  12.0  11.0   9.0   7.5

We swap 9.5 with 15.0 because it is the larger of the pair.

After we make the swap, we recursively call the heapIt function, starting at the swapped in position (position 5).

Position:     0     1     2     3     4     5     6     7     8     9    10    11
The Heap:  11.5   8.5  15.0  28.0  25.5   9.5  12.0  11.0   9.0   7.5    <>    <>

Here, we compare the value at position 5 versus 10 and 11, both of which are outside the bounds of our array.  So we do nothing. 

Position:     0     1     2     3     4     5     6     7     8     9    10    11
The Heap:  11.5   8.5  15.0  28.0  25.5   9.5  12.0  11.0   9.0   7.5    <>    <>

So, we pop back to where we were.  We need to move to the left and compare the values at position 1 versus 3 & 4. .  Here we have 8.5 versus the pair 28.0 and 25.5.  So, we make the swap for 28.0

Position:     0     1     2     3     4     5     6     7     8     9    10    11
The Heap:  11.5  28.0  15.0   8.5  25.5   9.5  12.0  11.0   9.0   7.5    <>    <>

After we make the swap, we recursively call the heapIt function, starting at the swapped in position (position 3), comparing to 7 & 8.

Position:     0     1     2     3     4     5     6     7     8     9    10    11
The Heap:  11.5  28.0  15.0   8.5  25.5   9.5  12.0  11.0   9.0   7.5    <>    <>

So, our first chore, is to take our floating point array of cars MPG's and transform it into a Heap using the rules above.

Here, 8.5 is trumped by 11.0, so it is (demoted) swapped out:

Position:     0     1     2     3     4     5     6     7     8     9    10    11
The Heap:  11.5  28.0  15.0  11.0  25.5   9.5  12.0   8.5   9.0   7.5    <>    <>

After we make the swap, we again recursively call the heapIt function, starting at the swapped in position (position 7), comparing to 15 & 16.  Obviously, there is nothing to do.

So, we again pop back to where we were.  We need to move to the left and compare the values at position 0 versus 1 & 2.  Here we have 11.5 versus the pair 28.0 and 15.0.  So, we make the swap for 28.0

Position:     0     1     2     3     4     5     6     7     8     9    10    11
The Heap:  11.5  28.0  15.0  11.0  25.5   9.5  12.0   8.5   9.0   7.5    <>    <>

Position:     0     1     2     3     4     5     6     7     8     9    10    11
The Heap:  28.0  11.5  15.0  11.0  25.5   9.5  12.0   8.5   9.0   7.5    <>    <>

After we make the swap, we again recursively call the heapIt function, starting at the swapped in position (position 1), comparing to 3 & 4.  This comparison yields a swap of position 1 for 4.

Position:     0     1     2     3     4     5     6     7     8     9    10    11
The Heap:  28.0  11.5  15.0  11.0  25.5   9.5  12.0   8.5   9.0   7.5    <>    <>

 This comparison yields a swap of position 1 for 4.

Position:     0     1     2     3     4     5     6     7     8     9    10    11
The Heap:  28.0  25.5  15.0  11.0  11.5   9.5  12.0   8.5   9.0   7.5    <>    <>

After we make the swap, we again recursively call the heapIt function, starting at the swapped in position (position 4), comparing to 9 & 10.  Because 11.5 is greater than the pair, we have nothing to do.

Our Heap is now done:

Position:     0     1     2     3     4     5     6     7     8     9    10    11
The Heap:  28.0  25.5  15.0  11.0  11.5   9.5  12.0   8.5   9.0   7.5    <>    <>

Now that we have the Heap, we are ready to process it.  This may seem strange, but this is how it works. 

1:  We know in a Heap the first element must be the biggest.  But, in our sorted array, we want it to be last.  What to do?

2:  Swap the first and last elements.  So we move the first to last and last to first.  We now have the last place settled.

In the example above we would have:

Position:     0     1     2     3     4     5     6     7     8     9
The Heap:   7.5  25.5  15.0  11.0  11.5   9.5  12.0   8.5   9.0  28.0

3.  We next shrink the working array by one; saving the last element because it is now set.

4.  We now want to check and see where the new value at position 0 fits.  We run our heapIt function to move it to the right.

I won't show all the steps, but we first compare element at zero with those at 1 and 2.  We swap 1 and 0, by our rules.  We then check the new position 1 against 3 & 4.  We swap 1 with 4.  We then check 4 against 9 & 10, but 9 & 10 don't exist (recall we shrunk our array by one), so we stop here:

Position:     0     1     2     3     4     5     6     7     8     9
The Heap:  25.5  11.5  15.0  11.0   7.5   9.5  12.0   8.5   9.0  28.0

5.  Once we stop, we again know that the first element is the largest.  We swap it with the tail (now position 8).  Again, we now want to check and see where the new value at position 0 fits.

6.  We continue this process until we precessed all our elements.

Here is the code:

-------------------------------------------------------

let rec ProcessHeap (rA : float array) (i : int) =
    match i with
    | 0 -> rA
    | _ -> let k = i-1
           let rB = swap rA 0 k
           let rA = heapIt rB 0 k
           printA i rA
           ProcessHeap rA k

let heapSort (rA : float array) (n : int) =
    // Initialize i = (n-1)/2
    // Less than or equal to midway
    // Have to account for zero-indexed array
    let i = (n-1)/2
    printA n rA
   
    //Build the Heap
    let rA = BuildTheHeap rA n i  
   
    for i = 0 to rA.Length-1 do
        printfn "The Heap: %10.3f" rA.[i]
    let rB = ProcessHeap rA n   
    rB

-------------------------------------------------------

Here we have our processHeap function that swaps the first and last elements and then reorganizes the Heap.  It also reduces the array bounds by one.  It then recursively calls itself, working element by element until the array size hits zero.

Our last function is just the generic heapSort function that takes in the array, does manages all the work and returns a newly reorganized array. 

I guess its time to stow this puppy into the Math Library and do some testing.  That will be my next post, but I'm tired and need a beer.

 

 

Print | posted on Tuesday, July 14, 2009 7:37 AM |

Powered by:
Powered By Subtext Powered By ASP.NET