Chorin's method is usually not a good choice for steady NS. Coupled solvers are highly recommended. Basically just add convection to one of the Stokes demos and use a nonlinear Newton solver because of the nonlinear convection. In short something like
v, q = TestFunctions(W)
up_ = Function(W)
u_, p_ = split(up_)
f = Constant((0, 0)) # Possible forcing term
nu = Constant(0.001)
F = inner(dot(grad(u_), u_), v)*dx + nu*inner(grad(u_), grad(v))*dx \
- inner(p_, div(v))*dx - inner(q, div(u_))*dx + inner(f, v)*dx
solve(F == 0, up_)