added some externally typed benchmarks as demos

7 files changed, 841 insertions, 0 deletions

diff --git a/Demos/benchmarks/nbody.pxd b/Demos/benchmarks/nbody.pxd
new file mode 100644
index 000000000..f58354fd4
--- /dev/null
+++ b/Demos/benchmarks/nbody.pxd
@@ -0,0 +1,18 @@
+
+cimport cython
+
+@cython.locals(x=Py_ssize_t)
+cdef combinations(list l)
+
+@cython.locals(x1=double, x2=double, y1=double, y2=double, z1=double, z2=double,
+               m1=double, m2=double, vx=double, vy=double, vz=double, i=long)
+cdef advance(double dt, long n, list bodies=*, list pairs=*)
+
+@cython.locals(x1=double, x2=double, y1=double, y2=double, z1=double, z2=double,
+               m1=double, m2=double, vx=double, vy=double, vz=double)
+cdef report_energy(list bodies=*, list pairs=*, double e=*)
+
+@cython.locals(vx=double, vy=double, vz=double, m=double)
+cdef offset_momentum(tuple ref, list bodies=*, double px=*, double py=*, double pz=*)
+
+cpdef test_nbody(long iterations)
diff --git a/Demos/benchmarks/nbody.py b/Demos/benchmarks/nbody.py
new file mode 100644
index 000000000..eee2a7477
--- /dev/null
+++ b/Demos/benchmarks/nbody.py
@@ -0,0 +1,150 @@
+#!/usr/bin/env python
+
+"""N-body benchmark from the Computer Language Benchmarks Game.
+
+This is intended to support Unladen Swallow's perf.py. Accordingly, it has been
+modified from the Shootout version:
+- Accept standard Unladen Swallow benchmark options.
+- Run report_energy()/advance() in a loop.
+- Reimplement itertools.combinations() to work with older Python versions.
+"""
+
+# Pulled from http://shootout.alioth.debian.org/u64q/benchmark.php?test=nbody&lang=python&id=4
+# Contributed by Kevin Carson.
+# Modified by Tupteq, Fredrik Johansson, and Daniel Nanz.
+
+__contact__ = "collinwinter@google.com (Collin Winter)"
+
+# Python imports
+import optparse
+import sys
+import time
+
+# Local imports
+import util
+
+def combinations(l):
+    """Pure-Python implementation of itertools.combinations(l, 2)."""
+    result = []
+    for x in xrange(len(l) - 1):
+        ls = l[x+1:]
+        for y in ls:
+            result.append((l[x],y))
+    return result
+
+
+PI = 3.14159265358979323
+SOLAR_MASS = 4 * PI * PI
+DAYS_PER_YEAR = 365.24
+
+BODIES = {
+    'sun': ([0.0, 0.0, 0.0], [0.0, 0.0, 0.0], SOLAR_MASS),
+
+    'jupiter': ([4.84143144246472090e+00,
+                 -1.16032004402742839e+00,
+                 -1.03622044471123109e-01],
+                [1.66007664274403694e-03 * DAYS_PER_YEAR,
+                 7.69901118419740425e-03 * DAYS_PER_YEAR,
+                 -6.90460016972063023e-05 * DAYS_PER_YEAR],
+                 9.54791938424326609e-04 * SOLAR_MASS),
+
+    'saturn': ([8.34336671824457987e+00,
+                4.12479856412430479e+00,
+                -4.03523417114321381e-01],
+               [-2.76742510726862411e-03 * DAYS_PER_YEAR,
+                4.99852801234917238e-03 * DAYS_PER_YEAR,
+                2.30417297573763929e-05 * DAYS_PER_YEAR],
+                2.85885980666130812e-04 * SOLAR_MASS),
+
+    'uranus': ([1.28943695621391310e+01,
+                -1.51111514016986312e+01,
+                -2.23307578892655734e-01],
+               [2.96460137564761618e-03 * DAYS_PER_YEAR,
+                2.37847173959480950e-03 * DAYS_PER_YEAR,
+                -2.96589568540237556e-05 * DAYS_PER_YEAR],
+                4.36624404335156298e-05 * SOLAR_MASS),
+
+    'neptune': ([1.53796971148509165e+01,
+                 -2.59193146099879641e+01,
+                 1.79258772950371181e-01],
+                [2.68067772490389322e-03 * DAYS_PER_YEAR,
+                 1.62824170038242295e-03 * DAYS_PER_YEAR,
+                 -9.51592254519715870e-05 * DAYS_PER_YEAR],
+                 5.15138902046611451e-05 * SOLAR_MASS) }
+
+
+SYSTEM = list(BODIES.values())
+PAIRS = combinations(SYSTEM)
+
+
+def advance(dt, n, bodies=SYSTEM, pairs=PAIRS):
+    for i in xrange(n):
+        for (([x1, y1, z1], v1, m1),
+             ([x2, y2, z2], v2, m2)) in pairs:
+            dx = x1 - x2
+            dy = y1 - y2
+            dz = z1 - z2
+            mag = dt * ((dx * dx + dy * dy + dz * dz) ** (-1.5))
+            b1m = m1 * mag
+            b2m = m2 * mag
+            v1[0] -= dx * b2m
+            v1[1] -= dy * b2m
+            v1[2] -= dz * b2m
+            v2[0] += dx * b1m
+            v2[1] += dy * b1m
+            v2[2] += dz * b1m
+        for (r, [vx, vy, vz], m) in bodies:
+            r[0] += dt * vx
+            r[1] += dt * vy
+            r[2] += dt * vz
+
+
+def report_energy(bodies=SYSTEM, pairs=PAIRS, e=0.0):
+    for (((x1, y1, z1), v1, m1),
+         ((x2, y2, z2), v2, m2)) in pairs:
+        dx = x1 - x2
+        dy = y1 - y2
+        dz = z1 - z2
+        e -= (m1 * m2) / ((dx * dx + dy * dy + dz * dz) ** 0.5)
+    for (r, [vx, vy, vz], m) in bodies:
+        e += m * (vx * vx + vy * vy + vz * vz) / 2.
+    return e
+
+
+def offset_momentum(ref, bodies=SYSTEM, px=0.0, py=0.0, pz=0.0):
+    for (r, [vx, vy, vz], m) in bodies:
+        px -= vx * m
+        py -= vy * m
+        pz -= vz * m
+    (r, v, m) = ref
+    v[0] = px / m
+    v[1] = py / m
+    v[2] = pz / m
+
+
+def test_nbody(iterations):
+    # Warm-up runs.
+    report_energy()
+    advance(0.01, 20000)
+    report_energy()
+
+    times = []
+    for _ in xrange(iterations):
+        t0 = time.time()
+        report_energy()
+        advance(0.01, 20000)
+        report_energy()
+        t1 = time.time()
+        times.append(t1 - t0)
+    return times
+
+
+if __name__ == '__main__':
+    parser = optparse.OptionParser(
+        usage="%prog [options]",
+        description=("Run the n-body benchmark."))
+    util.add_standard_options_to(parser)
+    options, args = parser.parse_args()
+
+    offset_momentum(BODIES['sun'])  # Set up global state
+    util.run_benchmark(options, options.num_runs, test_nbody)
diff --git a/Demos/benchmarks/richards.pxd b/Demos/benchmarks/richards.pxd
new file mode 100644
index 000000000..edc4df3f2
--- /dev/null
+++ b/Demos/benchmarks/richards.pxd
@@ -0,0 +1,91 @@
+cimport cython
+
+cdef class Packet:
+    cdef public object link
+    cdef public object ident
+    cdef public object kind
+    cdef public Py_ssize_t datum
+    cdef public list data
+
+    cpdef append_to(self,lst)
+
+cdef class TaskRec:
+    pass
+
+cdef class DeviceTaskRec(TaskRec):
+    cdef public object pending
+
+cdef class IdleTaskRec(TaskRec):
+    cdef public long control
+    cdef public Py_ssize_t count
+
+cdef class HandlerTaskRec(TaskRec):
+    cdef public object work_in   # = None
+    cdef public object device_in # = None
+
+    cpdef workInAdd(self,p)
+    cpdef deviceInAdd(self,p)
+
+cdef class WorkerTaskRec(TaskRec):
+    cdef public object destination # = I_HANDLERA
+    cdef public Py_ssize_t count
+
+cdef class TaskState:
+    cdef public bint packet_pending # = True
+    cdef public bint task_waiting   # = False
+    cdef public bint task_holding   # = False
+
+    cpdef packetPending(self)
+    cpdef waiting(self)
+    cpdef running(self)
+    cpdef waitingWithPacket(self)
+    cpdef bint isPacketPending(self)
+    cpdef bint isTaskWaiting(self)
+    cpdef bint isTaskHolding(self)
+    cpdef bint isTaskHoldingOrWaiting(self)
+    cpdef bint isWaitingWithPacket(self)
+
+cdef class TaskWorkArea:
+    cdef public list taskTab # = [None] * TASKTABSIZE
+
+    cdef public object taskList # = None
+
+    cdef public Py_ssize_t holdCount # = 0
+    cdef public Py_ssize_t qpktCount # = 0
+
+cdef class Task(TaskState):
+    cdef public Task link # = taskWorkArea.taskList
+    cdef public object ident # = i
+    cdef public object priority # = p
+    cdef public object input # = w
+    cdef public object handle # = r
+
+    cpdef addPacket(self,Packet p,old)
+    cpdef runTask(self)
+    cpdef waitTask(self)
+    cpdef hold(self)
+    cpdef release(self,i)
+    cpdef qpkt(self,Packet pkt)
+    cpdef findtcb(self,id)
+
+cdef class DeviceTask(Task):
+    @cython.locals(d=DeviceTaskRec)
+    cpdef fn(self,Packet pkt,r)
+
+cdef class HandlerTask(Task):
+    @cython.locals(h=HandlerTaskRec)
+    cpdef fn(self,Packet pkt,r)
+
+cdef class IdleTask(Task):
+    @cython.locals(i=IdleTaskRec)
+    cpdef fn(self,Packet pkt,r)
+
+cdef class WorkTask(Task):
+    @cython.locals(w=WorkerTaskRec)
+    cpdef fn(self,Packet pkt,r)
+
+@cython.locals(t=Task)
+cpdef schedule()
+
+cdef class Richards:
+    cpdef run(self, iterations)
diff --git a/Demos/benchmarks/richards.py b/Demos/benchmarks/richards.py
new file mode 100644
index 000000000..f718ad502
--- /dev/null
+++ b/Demos/benchmarks/richards.py
@@ -0,0 +1,440 @@
+# based on a Java version:
+#  Based on original version written in BCPL by Dr Martin Richards
+#  in 1981 at Cambridge University Computer Laboratory, England
+#  and a C++ version derived from a Smalltalk version written by
+#  L Peter Deutsch.
+#  Java version:  Copyright (C) 1995 Sun Microsystems, Inc.
+#  Translation from C++, Mario Wolczko
+#  Outer loop added by Alex Jacoby
+
+# Task IDs
+I_IDLE = 1
+I_WORK = 2
+I_HANDLERA = 3
+I_HANDLERB = 4
+I_DEVA = 5
+I_DEVB = 6
+
+# Packet types
+K_DEV = 1000
+K_WORK = 1001
+
+# Packet
+
+BUFSIZE = 4
+
+BUFSIZE_RANGE = range(BUFSIZE)
+
+class Packet(object):
+    def __init__(self,l,i,k):
+        self.link = l
+        self.ident = i
+        self.kind = k
+        self.datum = 0
+        self.data = [0] * BUFSIZE
+
+    def append_to(self,lst):
+        self.link = None
+        if lst is None:
+            return self
+        else:
+            p = lst
+            next = p.link
+            while next is not None:
+                p = next
+                next = p.link
+            p.link = self
+            return lst
+
+# Task Records
+
+class TaskRec(object):
+    pass
+
+class DeviceTaskRec(TaskRec):
+    def __init__(self):
+        self.pending = None
+
+class IdleTaskRec(TaskRec):
+    def __init__(self):
+        self.control = 1
+        self.count = 10000
+
+class HandlerTaskRec(TaskRec):
+    def __init__(self):
+        self.work_in = None
+        self.device_in = None
+
+    def workInAdd(self,p):
+        self.work_in = p.append_to(self.work_in)
+        return self.work_in
+
+    def deviceInAdd(self,p):
+        self.device_in = p.append_to(self.device_in)
+        return self.device_in
+
+class WorkerTaskRec(TaskRec):
+    def __init__(self):
+        self.destination = I_HANDLERA
+        self.count = 0
+# Task
+
+class TaskState(object):
+    def __init__(self):
+        self.packet_pending = True
+        self.task_waiting = False
+        self.task_holding = False
+
+    def packetPending(self):
+        self.packet_pending = True
+        self.task_waiting = False
+        self.task_holding = False
+        return self
+
+    def waiting(self):
+        self.packet_pending = False
+        self.task_waiting = True
+        self.task_holding = False
+        return self
+
+    def running(self):
+        self.packet_pending = False
+        self.task_waiting = False
+        self.task_holding = False
+        return self
+        
+    def waitingWithPacket(self):
+        self.packet_pending = True
+        self.task_waiting = True
+        self.task_holding = False
+        return self
+        
+    def isPacketPending(self):
+        return self.packet_pending
+
+    def isTaskWaiting(self):
+        return self.task_waiting
+
+    def isTaskHolding(self):
+        return self.task_holding
+
+    def isTaskHoldingOrWaiting(self):
+        return self.task_holding or (not self.packet_pending and self.task_waiting)
+
+    def isWaitingWithPacket(self):
+        return self.packet_pending and self.task_waiting and not self.task_holding
+
+
+
+
+
+tracing = False
+layout = 0
+
+def trace(a):
+    global layout
+    layout -= 1
+    if layout <= 0:
+        print
+        layout = 50
+    print a,
+
+
+TASKTABSIZE = 10
+
+class TaskWorkArea(object):
+    def __init__(self):
+        self.taskTab = [None] * TASKTABSIZE
+
+        self.taskList = None
+
+        self.holdCount = 0
+        self.qpktCount = 0
+
+taskWorkArea = TaskWorkArea()
+
+class Task(TaskState):
+
+
+    def __init__(self,i,p,w,initialState,r):
+        self.link = taskWorkArea.taskList
+        self.ident = i
+        self.priority = p
+        self.input = w
+
+        self.packet_pending = initialState.isPacketPending()
+        self.task_waiting = initialState.isTaskWaiting()
+        self.task_holding = initialState.isTaskHolding()
+
+        self.handle = r
+
+        taskWorkArea.taskList = self
+        taskWorkArea.taskTab[i] = self
+
+    def fn(self,pkt,r):
+        raise NotImplementedError
+
+
+    def addPacket(self,p,old):
+        if self.input is None:
+            self.input = p
+            self.packet_pending = True
+            if self.priority > old.priority:
+                return self
+        else:
+            p.append_to(self.input)
+        return old
+
+
+    def runTask(self):
+        if self.isWaitingWithPacket():
+            msg = self.input
+            self.input = msg.link
+            if self.input is None:
+                self.running()
+            else:
+                self.packetPending()
+        else:
+            msg = None
+
+        return self.fn(msg,self.handle)
+
+
+    def waitTask(self):
+        self.task_waiting = True
+        return self
+
+
+    def hold(self):
+        taskWorkArea.holdCount += 1
+        self.task_holding = True
+        return self.link
+
+
+    def release(self,i):
+        t = self.findtcb(i)
+        t.task_holding = False
+        if t.priority > self.priority:
+            return t
+        else:
+            return self
+
+
+    def qpkt(self,pkt):
+        t = self.findtcb(pkt.ident)
+        taskWorkArea.qpktCount += 1
+        pkt.link = None
+        pkt.ident = self.ident
+        return t.addPacket(pkt,self)
+
+
+    def findtcb(self,id):
+        t = taskWorkArea.taskTab[id]
+        if t is None:
+            raise Exception("Bad task id %d" % id)
+        return t
+            
+
+# DeviceTask
+
+
+class DeviceTask(Task):
+    def __init__(self,i,p,w,s,r):
+        Task.__init__(self,i,p,w,s,r)
+
+    def fn(self,pkt,r):
+        d = r
+        assert isinstance(d, DeviceTaskRec)
+        if pkt is None:
+            pkt = d.pending
+            if pkt is None:
+                return self.waitTask()
+            else:
+                d.pending = None
+                return self.qpkt(pkt)
+        else:
+            d.pending = pkt
+            if tracing: trace(pkt.datum)
+            return self.hold()
+
+
+
+class HandlerTask(Task):
+    def __init__(self,i,p,w,s,r):
+        Task.__init__(self,i,p,w,s,r)
+
+    def fn(self,pkt,r):
+        h = r
+        assert isinstance(h, HandlerTaskRec)
+        if pkt is not None:
+            if pkt.kind == K_WORK:
+                h.workInAdd(pkt)
+            else:
+                h.deviceInAdd(pkt)
+        work = h.work_in
+        if work is None:
+            return self.waitTask()
+        count = work.datum
+        if count >= BUFSIZE:
+            h.work_in = work.link
+            return self.qpkt(work)
+
+        dev = h.device_in
+        if dev is None:
+            return self.waitTask()
+
+        h.device_in = dev.link
+        dev.datum = work.data[count]
+        work.datum = count + 1
+        return self.qpkt(dev)
+
+# IdleTask
+
+
+class IdleTask(Task):
+    def __init__(self,i,p,w,s,r):
+        Task.__init__(self,i,0,None,s,r)
+
+    def fn(self,pkt,r):
+        i = r
+        assert isinstance(i, IdleTaskRec)
+        i.count -= 1
+        if i.count == 0:
+            return self.hold()
+        elif i.control & 1 == 0:
+            i.control /= 2
+            return self.release(I_DEVA)
+        else:
+            i.control = i.control/2 ^ 0xd008
+            return self.release(I_DEVB)
+            
+
+# WorkTask
+
+
+A = ord('A')
+
+class WorkTask(Task):
+    def __init__(self,i,p,w,s,r):
+        Task.__init__(self,i,p,w,s,r)
+
+    def fn(self,pkt,r):
+        w = r
+        assert isinstance(w, WorkerTaskRec)
+        if pkt is None:
+            return self.waitTask()
+
+        if w.destination == I_HANDLERA:
+            dest = I_HANDLERB
+        else:
+            dest = I_HANDLERA
+
+        w.destination = dest
+        pkt.ident = dest
+        pkt.datum = 0
+
+        for i in BUFSIZE_RANGE: # xrange(BUFSIZE)
+            w.count += 1
+            if w.count > 26:
+                w.count = 1
+            pkt.data[i] = A + w.count - 1
+
+        return self.qpkt(pkt)
+
+import time
+
+
+
+def schedule():
+    t = taskWorkArea.taskList
+    while t is not None:
+        pkt = None
+
+        if tracing:
+            print "tcb =",t.ident
+
+        if t.isTaskHoldingOrWaiting():
+            t = t.link
+        else:
+            if tracing: trace(chr(ord("0")+t.ident))
+            t = t.runTask()
+
+class Richards(object):
+
+    def run(self, iterations):
+        for i in xrange(iterations):
+            taskWorkArea.holdCount = 0
+            taskWorkArea.qpktCount = 0
+
+            IdleTask(I_IDLE, 1, 10000, TaskState().running(), IdleTaskRec())
+
+            wkq = Packet(None, 0, K_WORK)
+            wkq = Packet(wkq , 0, K_WORK)
+            WorkTask(I_WORK, 1000, wkq, TaskState().waitingWithPacket(), WorkerTaskRec())
+
+            wkq = Packet(None, I_DEVA, K_DEV)
+            wkq = Packet(wkq , I_DEVA, K_DEV)
+            wkq = Packet(wkq , I_DEVA, K_DEV)
+            HandlerTask(I_HANDLERA, 2000, wkq, TaskState().waitingWithPacket(), HandlerTaskRec())
+
+            wkq = Packet(None, I_DEVB, K_DEV)
+            wkq = Packet(wkq , I_DEVB, K_DEV)
+            wkq = Packet(wkq , I_DEVB, K_DEV)
+            HandlerTask(I_HANDLERB, 3000, wkq, TaskState().waitingWithPacket(), HandlerTaskRec())
+
+            wkq = None;
+            DeviceTask(I_DEVA, 4000, wkq, TaskState().waiting(), DeviceTaskRec());
+            DeviceTask(I_DEVB, 5000, wkq, TaskState().waiting(), DeviceTaskRec());
+            
+            schedule()
+
+            if taskWorkArea.holdCount == 9297 and taskWorkArea.qpktCount == 23246:
+                pass
+            else:
+                return False
+
+        return True
+
+def entry_point(iterations):
+    r = Richards()
+    startTime = time.time()
+    result = r.run(iterations)
+    endTime = time.time()
+    return result, startTime, endTime
+
+def main(entry_point = entry_point, iterations = 10):
+    print "Richards benchmark (Python) starting... [%r]" % entry_point
+    result, startTime, endTime = entry_point(iterations)
+    if not result:
+        print "Incorrect results!"
+        return -1
+    print "finished."
+    total_s = endTime - startTime
+    print "Total time for %d iterations: %.2f secs" %(iterations,total_s)
+    print "Average time per iteration: %.2f ms" %(total_s*1000/iterations)
+    return 42
+
+try:
+    import sys
+    if '-nojit' in sys.argv:
+        sys.argv.remove('-nojit')
+        raise ImportError
+    import pypyjit
+except ImportError:
+    pass
+else:
+    import types
+    for item in globals().values():
+        if isinstance(item, types.FunctionType):
+            pypyjit.enable(item.func_code)
+        elif isinstance(item, type):
+            for it in item.__dict__.values():
+                if isinstance(it, types.FunctionType):
+                    pypyjit.enable(it.func_code)
+    
+if __name__ == '__main__':
+    import sys
+    if len(sys.argv) >= 2:
+        main(iterations = int(sys.argv[1]))
+    else:
+        main()
diff --git a/Demos/benchmarks/spectralnorm.pxd b/Demos/benchmarks/spectralnorm.pxd
new file mode 100644
index 000000000..827826d5d
--- /dev/null
+++ b/Demos/benchmarks/spectralnorm.pxd
@@ -0,0 +1,18 @@
+
+cimport cython
+
+cdef inline double eval_A(double i, double j)
+
+@cython.locals(i=long)
+cdef list eval_A_times_u(list u)
+
+@cython.locals(i=long)
+cdef list eval_At_times_u(list u)
+
+cdef list eval_AtA_times_u(list u)
+
+@cython.locals(j=long, u_j=double, partial_sum=double)
+cdef double part_A_times_u(double i, list u)
+
+@cython.locals(j=long, u_j=double, partial_sum=double)
+cdef double part_At_times_u(double i, list u)
diff --git a/Demos/benchmarks/spectralnorm.py b/Demos/benchmarks/spectralnorm.py
new file mode 100644
index 000000000..abbb32b69
--- /dev/null
+++ b/Demos/benchmarks/spectralnorm.py
@@ -0,0 +1,69 @@
+# -*- coding: utf-8 -*-
+# The Computer Language Benchmarks Game
+# http://shootout.alioth.debian.org/
+# Contributed by Sebastien Loisel
+# Fixed by Isaac Gouy
+# Sped up by Josh Goldfoot
+# Dirtily sped up by Simon Descarpentries
+# Concurrency by Jason Stitt
+
+from math            import sqrt
+from itertools       import izip
+from time import time
+import util
+import itertools
+import optparse
+
+def eval_A (i, j):
+    return 1.0 / ((i + j) * (i + j + 1) / 2 + i + 1)
+
+def eval_A_times_u (u):
+    return [ part_A_times_u(i,u) for i in xrange(len(u)) ]
+
+def eval_At_times_u (u):
+    return [ part_At_times_u(i,u) for i in xrange(len(u)) ]
+
+def eval_AtA_times_u (u):
+    return eval_At_times_u (eval_A_times_u (u))
+
+def part_A_times_u(i, u):
+    partial_sum = 0
+    for j, u_j in enumerate(u):
+        partial_sum += eval_A (i, j) * u_j
+    return partial_sum
+
+def part_At_times_u(i, u):
+    partial_sum = 0
+    for j, u_j in enumerate(u):
+        partial_sum += eval_A (j, i) * u_j
+    return partial_sum
+
+DEFAULT_N = 130
+
+def main(n):
+    times = []
+    for i in range(n):
+        t0 = time()
+        u = [1] * DEFAULT_N
+
+        for dummy in xrange (10):
+            v = eval_AtA_times_u (u)
+            u = eval_AtA_times_u (v)
+
+        vBv = vv = 0
+
+        for ue, ve in izip (u, v):
+            vBv += ue * ve
+            vv  += ve * ve
+        tk = time()
+        times.append(tk - t0)
+    return times
+    
+if __name__ == "__main__":
+    parser = optparse.OptionParser(
+        usage="%prog [options]",
+        description="Test the performance of the spectralnorm benchmark")
+    util.add_standard_options_to(parser)
+    options, args = parser.parse_args()
+
+    util.run_benchmark(options, options.num_runs, main)
diff --git a/Demos/benchmarks/util.py b/Demos/benchmarks/util.py
new file mode 100644
index 000000000..95495e713
--- /dev/null
+++ b/Demos/benchmarks/util.py
@@ -0,0 +1,55 @@
+#!/usr/bin/env python
+
+"""Utility code for benchmark scripts."""
+
+__author__ = "collinwinter@google.com (Collin Winter)"
+
+import math
+import operator
+
+try:
+    reduce
+except NameError:
+    from functools import reduce
+
+def run_benchmark(options, num_runs, bench_func, *args):
+    """Run the given benchmark, print results to stdout.
+
+    Args:
+        options: optparse.Values instance.
+        num_runs: number of times to run the benchmark
+        bench_func: benchmark function. `num_runs, *args` will be passed to this
+            function. This should return a list of floats (benchmark execution
+            times).
+    """
+    if options.profile:
+        import cProfile
+        prof = cProfile.Profile()
+        prof.runcall(bench_func, num_runs, *args)
+        prof.print_stats(sort=options.profile_sort)
+    else:
+        data = bench_func(num_runs, *args)
+        if options.take_geo_mean:
+            product = reduce(operator.mul, data, 1)
+            print(math.pow(product, 1.0 / len(data)))
+        else:
+            for x in data:
+                print(x)
+
+
+def add_standard_options_to(parser):
+    """Add a bunch of common command-line flags to an existing OptionParser.
+
+    This function operates on `parser` in-place.
+
+    Args:
+        parser: optparse.OptionParser instance.
+    """
+    parser.add_option("-n", action="store", type="int", default=100,
+                      dest="num_runs", help="Number of times to run the test.")
+    parser.add_option("--profile", action="store_true",
+                      help="Run the benchmark through cProfile.")
+    parser.add_option("--profile_sort", action="store", type="str",
+                      default="time", help="Column to sort cProfile output by.")
+    parser.add_option("--take_geo_mean", action="store_true",
+                      help="Return the geo mean, rather than individual data.")