  We present exact methods for the numerical integration of the Wannier-Stark
system in a many-body scenario including two Bloch bands. Our ab initio
approaches allow for the treatment of a few-body problem with bosonic
statistics and strong interparticle interaction. The numerical implementation
is based on the Lanczos algorithm for the diagonalization of large, but sparse
symmetric Floquet matrices. We analyze the scheme efficiency in terms of the
computational time, which is shown to scale polynomially with the size of the
system. The numerically computed eigensystem is applied to the analysis of the
Floquet Hamiltonian describing our problem. We show that this allows, for
instance, for the efficient detection and characterization of avoided crossings
and their statistical analysis. We finally compare the efficiency of our
Lanczos diagonalization for computing the temporal evolution of our many-body
system with an explicit fourth order Runge-Kutta integration. Both
implementations heavily exploit efficient matrix-vector multiplication schemes.
Our results should permit an extrapolation of the applicability of exact
methods to increasing sizes of generic many-body quantum problems with bosonic
statistics.
