  Kinetic theories constitute one of the most promising tools to decipher the
characteristic spatio-temporal dynamics in systems of actively propelled
particles. In this context, the Boltzmann equation plays a pivotal role, since
it provides a natural translation between a particle-level description of the
system's dynamics and the corresponding hydrodynamic fields. Yet, the intricate
mathematical structure of the Boltzmann equation substantially limits the
progress toward a full understanding of this equation by solely analytical
means. Here, we propose a general framework to numerically solve the Boltzmann
equation for self-propelled particle systems in two spatial dimensions and with
arbitrary boundary conditions. We discuss potential applications of this
numerical framework to active matter systems, and use the algorithm to give a
detailed analysis to a model system of self-propelled particles with polar
interactions. In accordance with previous studies, we find that spatially
homogeneous isotropic and broken symmetry states populate two distinct regions
in parameter space, which are separated by a narrow region of spatially
inhomogeneous, density-segregated moving patterns. We find clear evidence that
these three regions in parameter space are connected by first order phase
transitions, and that the transition between the spatially homogeneous
isotropic and polar ordered phases bears striking similarities to liquid-gas
phase transitions in equilibrium systems. Within the density segregated
parameter regime, we find a novel stable limit-cycle solution of the Boltzmann
equation, which consists of parallel lanes of polar clusters moving in opposite
directions, so as to render the overall symmetry of the system's ordered state
nematic, despite purely polar interactions on the level of single particles.
