Longest Increasing Subsequence, O(n log(n))

#!/usr/bin/env python
# @ hieuza [x] gmail.com

# longest increasing sequence
# O(n log n) algorithm
# http://en.wikipedia.org/wiki/Longest_increasing_subsequence

def simple_search(s, i, m, p, L):
    # simple linear search
    for j in xrange(L, 0, -1):
        if s[m[j]] < s[i]:
            return j
    return 0


# binary search
# search for *max* j in lo..hi, s.t., s[m[j]] < s[i]
def binary_search(s, i, m, p, lo, hi):
    if s[m[hi]] < s[i]:
        return hi

    if lo > hi:
        return 0 # not found

    mid = lo + (hi - lo)/2
    if s[m[mid]] < s[i]:
        return binary_search(s, i, m, p, mid, hi - 1)
    else:
        return binary_search(s, i, m, p, lo, mid - 1)


def search(s, i, m, p, L):
    return binary_search(s, i, m, p, 1, L)


def extract_lis(s, m, p, L):
    lis = []
    x = m[L]
    for _ in xrange(L, 0, -1):
        lis = [s[x]] + lis
        x = p[x]

    return lis


def lis(s):
    n = len(s)

    # for an increasing sequence of length j,
    # m[j] = k is the position of the smallest value s[k]
    # that is end of the sequence of length j.
    # so s[m[1]], s[m[2]], ... s[m[L]] is the longest sequence.
    m = [-1]*(n + 1)

    # p[i] is the position of the predecessor of s[i] in an increasing sequence.
    # need to use p[i] in binary search
    p = [0]*n

    # current max length
    L = 0

    for i, si in enumerate(s):
        # search for maximum j, s.t., s[m[j]] < s[i]
        # i.e., j is the max length of an i.s that has end value less than s[i]
        j = search(s, i, m, p, L)

        # parrent of i
        p[i] = m[j]

        if j == L or s[m[j + 1]] > s[i]:
            # update the end of sequence of length j + 1 to be i
            # it can happen when j == L, i.e., x[i] extends the current longest i.s
            # OR, x[i] is a new end of an existing i.s. of length j + 1
            m[j + 1] = i

            L = max(L, j + 1)

    # only print small list
    if n <= 30:
        print 's:', s
        print 'm:', m
        print 'p:', p
        lis = extract_lis(s, m, p, L)
        print 'lis:', lis

    print 'max LIS length:', L


def test(n, vmax):
    import random
    # s = [1, 2, 1]
    s = [random.randint(1, vmax) for _ in xrange(n)]
    lis(s)


if __name__ == '__main__':
    import sys

    n = 10
    vmax = 20

    print 'usage: [n] [vmax]'

    if len(sys.argv) > 1:
        n = int(sys.argv[1])

    if len(sys.argv) > 2:
        vmax = int(sys.argv[2])

    test(n, vmax)