Using 'yield' to simulate a Markov chain

One thing that I really like about Sage is that it uses Python as its underlying language. This means that we get “for free” many nice features of Python. One of these features that I particularly like is the yield keyword. Here is a small example:

def foo():
    i = 0
    while True:
        yield i
        i = i + 1

We can use the foo function as a generator:

sage: g = foo()
sage: print g.next()
0
sage: print g.next()
1
sage: print g.next()
2

In other words, the yield keyword acts as a way to a Python function into a generator. The execution of foo is paused until the next call to g.next(). If we reach the end of the function, the StopIteration exception is raised.

The yield keyword makes it pretty easy to write the skeleton for a Markov chain simulator, using the following basic form:

def markov_chain():
    state = initial_state()
    while True:
        yield state

        state = new_state(state)

        if some_condition: return

For a real example, see sage.combinat.matrices.latin.LatinSquare_generator. The Markov chain itself was given by Jacobson and Matthews in “Generating uniformly distributed random Latin squares,” Journal of Combinatorial Designs, vol 4, 1996, no 6, pp 405–437.