We present the theoretical foundations and the implementation details of a density-functional approach for coupled photons, electrons, and effective nuclei in non-relativistic quantum electrodynamics. Starting point of the formalism is a generalization of the Pauli-Fierz field theory for which we establish a one-to-one correspondence between external fields and internal variables. Based on this correspondence, we introduce a Kohn-Sham construction which provides a computationally feasible approach for ab-initio light-matter interactions. In the mean-field limit for the effective nuclei the formalism reduces to coupled Ehrenfest-Maxwell-Pauli-Kohn-Sham equations. We present an implementation of the approach in the real-space real-time code Octopus. For the implementation we use the Riemann-Silberstein formulation of classical electrodynamics and rewrite Maxwell's equations in Schrödinger form. This allows us to use existing time-evolution algorithms developed for quantum-mechanical systems also for Maxwell's equations. We introduce a predictorcorrector scheme and show how to couple the Riemann-Silberstein time-evolution of the electromagnetic fields self-consistently to the time-evolution of the electrons and nuclei. Furthermore, the Riemann-Silberstein approach allows to seamlessly combine macroscopic dielectric media with a microscopic coupling to matter currents. For an efficient absorption of outgoing electromagnetic waves, we present a perfectly matched layer for the Riemann-Silberstein vector. We introduce the concept of electromagnetic detectors, which allow to measure outgoing radiation in the far field and provide a direct way to record various spectroscopies. We present a multi-scale approach in space and time which allows to deal with the different length-scales of light and matter for a multitude of applications. We apply the approach to laser-induced plasmon excitation in a nanoplasmonic dimer system. We find that the self-consistent coupling of light and matter leads to significant local field effects which can not be captured with the conventional light-matter forward coupling. For our nanoplasmonic example we show that the self-consistent foward-backward coupling leads to changes in observables which are larger than the difference between local density and gradient corrected approximations for the exchange correlation functional. In addition, in our example we observe harmonic generation which appears only beyond the dipole approximation and can be directly observed in the outgoing electromagnetic waves on the simulation grid. The self-consistent coupling of the electromagnetic fields to the ion motion reveals significant energy transfer from the electromagnetic fields to matter on the scale of a few tens of femtoseconds. Overall, our approach is ideally suited for applications in nano-optics, nano-plasmonics, (photo) electrocatalysis, light-matter coupling in 2D materials, cases where laser pulses carry orbital angular momentum, or light-tailored chemical reactions in optical cavities to name but ...