Mathematica to Python

How can this Mathematica code be ported to Python? I do not know the Mathematica syntax and am having a hard time understanding how this is described in a more traditional language.

Source (pg 5): http://subjoin.net/misc/m496pres1.nb.pdf

This cannot be ported to Python directly as the definition a[j] uses the Symbolic Arithmetic feature of Mathematica.

a[j] is basically the coefficient of x^j in the series expansion of that rational function inside Apart.

Assume you have a[j], then f[n] is easy. A Block in Mathematica basically introduces a scope for variables. The first list initializes the variable, and the rest is the execution of the code. So

from __future__ import division
def f(n):
  v = n // 5
  q = v // 20
  r = v % 20
  return sum(binomial(q+5-j, 5) * a[r+20*j] for j in range(5))

(binomial is the Binomial coefficient.)

Using the proposed solutions from the previous answers I found that sympy sadly doesn't compute the apart() of the rational immediatly. It somehow gets confused. Moreover, the python list of coefficients returned by *Poly.all_coeffs()* has a different semantics than a Mathmatica list. Hence the try-except-clause in the definition of a().

The following code does work and the output, for some tested values, concurs with the answers given by the Mathematica formula in Mathematica 7:

from __future__ import division
from sympy import expand, Poly, binomial, apart
from sympy.abc import x

A = Poly(apart(expand(((1-x**20)**5)) / expand((((1-x)**2)*(1-x**2)*(1-x**5)*(1-x**10))))).all_coeffs()

def a(n):
    try:
        return A[n]
    except IndexError:
        return 0

def f(n):
    v = n // 5
    q = v // 20
    r = v % 20
    return sum(a[r+20*j]* binomial(q+5-j, 5) for j in range(5))

print map(f, [100, 50, 1000, 150])

The symbolics can be done with sympy. Combined with KennyTM's answer, something like this might be what you want:

from __future__ import division
from sympy import Symbol, apart, binomial

x = Symbol('x')
poly = (1-x**20)**5 / ((1-x)**2 * (1-x**2) * (1-x**5) * (1-x**10))
poly2 = apart(poly,x)

def a(j):
    return poly2.coeff(x**j)

def f(n):
    v = n // 5
    q = v // 20
    r = v % 20
    return sum(binomial(q+5-j, 5)*a(r+20*j) for j in range(5))

Although I have to admit that f(n) does not work (I'm not very good at Python).