Algorithm Implementation/Sorting/Patience sort

Given a generic STL-style input iterator to a sequence, we shall call its elements' value type card_type, and assume that we have some container of cards PileType (an example could be vector<card_type>), and a container of piles TableType (which, again, could be just a vector of piles).

In this case, one can define the patience sorting as follows:

    /**
     * Patience sort algorithm as per Aldous and Diaconis,
     * using binary search to select the leftmost heap for a new card.
     * Compares cards using operator<
     * 
     * TableType requirements
     *  -# whatever lower_bound needs
     *  -# begin(), end(): ForwardIterator of piles
     * 
     * PileType requirements
     *  -# default constructible
     *  -# constructible from a single card const reference
     *      typedef typename PileType::value_type card_type;
     *  -# push_back(const card_type&)
     *  -# assignable, lightweight on a singleton pile
     *  -# whatever lower_bound needs
     *  -# back() implemented or an appropriate pile_less specialization
     */
    template < typename PileType, typename TableType, typename InputIterator > 
    void patience_sort(TableType& table, InputIterator first, InputIterator last) {
        typedef PileType pile_type;
        typedef typename PileType::value_type card_type;
        typedef TableType table_type;
        typedef typename table_type::iterator iterator;
  
        while (first != last) {
            pile_type new_pile(*first); // read *first only once! (input iterator)
            // find the leftmost heap with the top card greater 
            // than *first , to push *first on it
            iterator pile_p = std::lower_bound(
                                table.begin(), table.end(), new_pile, 
                                pile_less<pile_type>() );
            if (pile_p != table.end()) {
                pile_p->push_back(new_pile.back());
            }
            // ...or allocate a new heap on the right with *first
            else {
                table.push_back(new_pile);
            }
            first++;
        }
    }

For the std::lower_bound algorithm to work, a tiny utility comparator functor is needed:

  /**
   * Comparator for two piles according to their top card values.
   * Compares cards using operator<
   */
  template < typename pile_type >
  struct pile_less : std::binary_function<pile_type, pile_type, bool>
  {
      bool operator()(const pile_type& x, const pile_type& y) const
      {
          return x.back() < y.back();
      }
  };
/**
 * When the algorithm is run with the sole purpose to discover the length of the longest increasing subsequence, 
 * one doesn't need to store any cards in the piles beneath the topmost one, since only the topmost card's value
 * is ever used by the algorithm above, and the length of the longest increasing subsequence is equal to the number of piles on the
 * table,
 * i.e., ''table.size()'':
 */
    /** 
     * A class representing a pile that is sticky, i.e., once a card
     * is put onto its top, it's glued forever. This way, only
     * the top card is actually remembered (O(1) memory per heap).
     * Used by longest_increasing_subsequence().
     *
     * Unfortunately, such heap is impossible to use in more advanced
     * algorithms, such as the Bespamyatnykh and Segal one, for the
     * actual longest increasing subsequence recovery.
     */
    template< typename Card >
    struct StickyPile { // cards stick together when put here :-)
        typedef Card value_type;
        typedef const value_type& const_reference;
        StickyPile(const Card& _card) : card(_card) {}
        void push_back(const Card& _card) { card = _card; }
        const_reference back() const { return card; }
    private:
        Card card;
    };
  
    /** Find the length of the longest increasing subsequence
     * given by the given iterator range, using patience_sort().
     * Compares cards using operator<
     */
    template< typename InputIterator >
    size_t longest_increasing_subsequence( 
            InputIterator first, InputIterator last)
    {
        typedef typename std::iterator_traits<InputIterator>::value_type Card;
        typedef StickyPile<Card> pile_type;
        typedef std::vector< pile_type > table_type;
        table_type table;
        patience_sort<pile_type>(table, first, last);
        return table.size();
    }

/* The above GPL code is based on excerpts from [[http://www.tarunz.org/~vassilii/pub/c++/patience_sort.hpp patience_sort.hpp]],
 which also provides versions of the same functionality parametrized by any given comparator (rather than using operator< ). */

Java

edit
import java.util.*;
public class PatienceSort
{
    public static <E extends Comparable<? super E>> void sort (E[] n)
    {
        List<Pile<E>> piles = new ArrayList<Pile<E>>();
        // sort into piles
        for (E x : n)
        {
            Pile<E> newPile = new Pile<E>();
            newPile.push(x);
            int i = Collections.binarySearch(piles, newPile);
            if (i < 0) i = ~i;
            if (i != piles.size())
                piles.get(i).push(x);
            else
                piles.add(newPile);
        }
        System.out.println("longest increasing subsequence has length = " + piles.size());
        
        // priority queue allows us to retrieve least pile efficiently
        PriorityQueue<Pile<E>> heap = new PriorityQueue<Pile<E>>(piles);
        for (int c = 0; c < n.length; c++)
        {
            Pile<E> smallPile = heap.poll();
            n[c] = smallPile.pop();
            if (!smallPile.isEmpty())
                heap.offer(smallPile);
        }
        assert(heap.isEmpty());
    }
    
    private static class Pile<E extends Comparable<? super E>> extends Stack<E> implements Comparable<Pile<E>>
    {
        public int compareTo(Pile<E> y) { return peek().compareTo(y.peek()); }
    }
}

Python 3

edit
import bisect
import heapq


def sort(seq):
    piles = []
    for x in seq:
        new_pile = [x]
        i = bisect.bisect_left(piles, new_pile)
        if i != len(piles):
            piles[i].insert(0, x)
        else:
            piles.append(new_pile)
    print("longest increasing subsequence has length =", len(piles))

    # priority queue allows us to retrieve least pile efficiently
    for i in range(len(seq)):
        small_pile = piles[0]
        seq[i] = small_pile.pop(0)
        if small_pile:
            heapq.heapreplace(piles, small_pile)
        else:
            heapq.heappop(piles)
    assert not piles


foo = [4, 65, 2, 4, -31, 0, 99, 1, 83, 782, 1]
sort(foo)
print(foo)

Clojure

edit

Implementation using the Patience Sort approach. The elements (newelem) put on a pile combine the "card" with a reference to the top of the previous stack, as per the algorithm. The combination is done using cons, so what gets put on a pile is a list—a descending subsequence.[1]

(defn place [piles card]
  (let [[les gts] (->> piles (split-with #(<= (ffirst %) card)))
        newelem (cons card (->> les last first))
        modpile (cons newelem (first gts))]
    (concat les (cons modpile (rest gts)))))
 
(defn a-longest [cards]
  (let [piles (reduce place '() cards)]
    (->> piles last first reverse)))
 
(println (a-longest [3 2 6 4 5 1]))
(println (a-longest [0 8 4 12 2 10 6 14 1 9 5 13 3 11 7 15]))

Output:

(2 4 5)
(0 2 6 9 11 15)
package main

import (
    "fmt"
    "container/heap"
)

type PileHeap [][]int

func (h PileHeap) Len() int      { return len(h) }
func (h PileHeap) Swap(i, j int) { h[i], h[j] = h[j], h[i] }
func (h PileHeap) Less(i, j int) bool {
    return h[i][len(h[i])-1] < h[j][len(h[j])-1]
}

func (h *PileHeap) Push(x interface{}) {
    *h = append(*h, x.([]int))
}

func (h *PileHeap) Pop() interface{} {
    old := *h
    n := len(old)
    x := old[n-1]
    *h = old[0 : n-1]
    return x
}

/*
bisectPilesRight uses binary search to returns the index where to insert card x,
assuming piles is already sorted according to the value of the top card
in each pile

The return value i is such that it's the largest i
for which the top card in piles[i] >= x and
return i == len(piles) if no such pile can be found
*/
func bisectPilesRight(piles [][]int , x int) int {
    lo, hi := 0, len(piles)
    for lo < hi {
        // invariant: x maybe between a[lo]...a[hi-1]
        mid := lo + (hi-lo)/2  // don't use (lo+hi)/2 to avoid overflow
        // Note that since (hi-lo)/2 >= 0, lo <= mid < hi
        pile := piles[mid]
        if x < pile[len(pile) - 1] { // compare x to top of pile
            hi = mid   // x may be between a[lo]...a[mid-1]
        } else {
            lo = mid+1 // x may be between a[mid+1]...a[hi]
        }
        // The new range is either lo...mid or mid+1...hi and
        // because lo<=mid<hi, the new range is always smaller than lo..hi
    }
    return lo
}

func PatienceSort(a []int) []int {
    piles := make([][]int, 0, 10)  // each pile will be a slice. 
    for _, x := range(a) {
        i := bisectPilesRight(piles, x)
        if i < len(piles) {
            piles[i] = append(piles[i], x)
        } else {
            piles = append(piles, []int{x}) // make a new pile
        }
        // fmt.Println(piles)
    }

    h := PileHeap(piles) // Use piles as a heap
    // heap.Init(&h) is not need because piles are already sorted by top card
    n := len(a)
    sorted := make([]int, n)
    for i := 0; i < n; i++ {
        pile := heap.Pop(&h).([]int)
        top := len(pile) - 1
        sorted[i] = pile[top]
        if top > 0 {
            // Put pile minus the top card back in heap if it is not empty
            heap.Push(&h, pile[:top])
        }
    }
    return sorted
}

func main() {
    a := []int{2,6,3,1,5,9,2}
    fmt.Print(patienceSort(a))
}

Phix

edit
function patience_sort(sequence s)
    -- create list of sorted lists
    sequence piles = {}
    for i=1 to length(s) do
        object n = s[i]
        for p=1 to length(piles)+1 do
            if p>length(piles) then
                piles = append(piles,{n})
            elsif n>=piles[p][$] then
                piles[p] = append(piles[p],n)
                exit
            end if
        end for
    end for
    -- merge sort the piles
    sequence res = ""
    while length(piles) do
        integer idx = smallest(piles,return_index:=true)
        res = append(res,piles[idx][1])
        if length(piles[idx])=1 then
            piles[idx..idx] = {}
        else
            piles[idx] = piles[idx][2..$]
        end if
    end while
    return res
end function


References

edit