The Maxwell equations in the unbounded three dimensional space are coupled to the Landau-Lifshitz-Gilbert equation on a (not necessarily convex) bounded domain. A weak formulation of the whole coupled system is derived based on the boundary integral formulation of the exterior Maxwell equations. We show existence of a weak solution and uniqueness of the Maxwell part of the weak solution. A numerical algorithm is proposed based on finite elements and boundary elements as spatial discretisation and using the backward Euler method and convolution quadratures for the interior domain and the boundary, respectively. Well-posedness and convergence of the numerical algorithm are shown, under minimal assumptions on the regularity of solutions. Numerical experiments illustrate and expand on the theoretical results.