We introduce a boundary integral method for two-dimensional quantum billiards subjected to a constant magnetic field. It allows to calculate spectra and wave functions, in particular at strong fields and semiclassical values of the magnetic length. The method is presented for interior and exterior problems with general boundary conditions. We explain why the magnetic analogues of the field-free single and double layer equations exhibit an infinity of spurious solutions and how these can be eliminated at the expense of dealing with (hyper-)singular operators. The high efficiency of the method is demonstrated by numerical calculations in the extreme semiclassical regime.