In this paper we present some computational techniques based on the class of preconditioned Krylov subspace methods that enable us to carry out large-scale, big data simulations of Computational Electromagnetics applications modeled using integral equations. This analysis requires the solution of large linear systems that cannot be afforded by conventional direct methods (based on variants of the Gaussian elimination algorithm) due to their high memory costs. We show that, thanks to the development of efficient Krylov methods and suitable preconditioning techniques, nowadays the solution of realistic electromagnetic problems that involve tens of million (and sometimes even more) unknowns, has become feasible. However, the choice of the best class of methods for the selected computer hardware and the given geometry remains an open problem that requires further analysis.