Although losses in metal are viewed as a serious drawback in many plasmonics experiments, thermoplasmonics is the field of physics that tries to take advantage of the latter. Indeed, the strong field enhancement obtained in nanometallic structures lead to a localized temperature increase in its vicinity, leading to interesting photothermal effects. Therefore, metallic nanoparticles may be used as heat sources that can be easily integrated in various environments. This is especially appealing in the field of nanomedecine and can for example be used for diagnosis purposes or nanosurgery to cite but just a few. Due to the various scales and phenomena that come into play, accurate numerical modeling is challenging. Laser illumination first excites a plasmon oscillation (reaction of the electrons of the metal) that relaxes to a thermal equilibrium and in turn excites the metal lattice (phonons). The latter is then responsible for heating the surroundings. A relevant modeling approach thus consists in describing the electron-phonon coupling through the evolution of their respective temperature. Maxwell's equations are then coupled to a set of coupled nonlinear hyperbolic equations for the temperatures of respectively electrons, phonons and environment. The nonlinearities and the different time scales at which each thermalization occurs make the numerical approximation of these equations quite challenging. We develop a suitable numerical framework based on high order discontinuous Galerkin methods for studying thermoplasmonics in the time-domain setting.