We utilize a Fourier transformation-based representation of Maxwell's equations to develop physics-constrained neural networks (PCNN) for electrodynamics without gauge ambiguity, which we label the Fourier-Helmholtz-Maxwell Neural Operator (FoHM-NO) method. In this approach, both of Gauss's laws and Faraday's law are built in as hard constraints, as well as the longitudinal component of Ampère-Maxwell in Fourier space, assuming the continuity equation. An encoder-decoder network acts as a solution operator for the transverse components of the Fourier transformed vector potential, A⟂,z(k,t), whose two degrees of freedom are used to predict the electromagnetic fields. This method was tested on two electron beam simulations. Among the models investigated, it was found that a U-Net architecture exhibited the best performance as it trained quicker, was more accurate and generalized better than the other architectures examined. We demonstrate that our approach is useful for solving Maxwell's equations for the electromagnetic fields generated by intense relativistic charged particle beams and that it generalizes well to unseen test data, while being orders of magnitude quicker than conventional simulations. We show that the model can be re-trained to make highly accurate predictions in as few as 20 epochs on a previously unseen data set.