We numerically investigate the dynamics near black hole formation of solutions to the Einstein-Vlasov system in axisymmetry. A symmetry reduction of the Einstein-Vlasov equations in the (2+1)+1 formalism is performed, including rotation. The numerical implementation is based on a constrained (mixed hyperbolic-elliptic) evolution scheme. The collisionless matter is treated using a particle-in-cell method. We tune smooth one-parameter families of initial data to the threshold of black hole formation and observe type I critical behaviour, in particular power-law scaling of the lifetime of the near-critical solution. The qualitative behaviour close to the critical point is found to depend on the sign of the binding energy. We prove that complete dispersal of the solution implies that it has nonpositive binding energy. This is joint work with Ellery Ames and Håkan Andréasson.