We present a family of unconditionally stable algorithms, based on the Suzuki product-formula approach, that solve the time-dependent Maxwell equations in systems with spatially varying permittivity and permeability. Salient features of these algorithms are discussed. As an illustration we compute, as a function of cluster size, the spectrum of electromagnetic modes in a cluster of photonic bandgap material