A comprehensive study of the atmospheric neutrino flux in the energy region from sub-GeV up to several TeV using the Super-Kamiokande (SK) water Cherenkov detector is presented in this paper. The energy and azimuthal spectra, and variation over time, of the atmospheric νe+νe and νμ+νμ fluxes are measured. The energy spectra are obtained using an iterative unfolding method by combining various event topologies with differing energy responses. The azimuthal spectra depending on energy and zenith angle, and their modulation by geomagnetic effects, are also studied. A predicted east-west asymmetry is observed in both the νe and νμ samples at 8.0σ and 6.0σ significance, respectively, and an indication that the asymmetry dipole angle changes depending on the zenith angle was seen at the 2.2σ level. The measured energy and azimuthal spectra are consistent with the current flux models within the estimated systematic uncertainties. A study of the long-term correlation between the atmospheric neutrino flux and the solar magnetic activity cycle is performed, and a weak preference for a correlation was seen at the 1.1σ level, using SK-I-SK-IV data spanning a 20-year period. For several particularly strong solar activity periods, corresponding to Forbush decrease events, no theoretical prediction is available but a deviation below the typical neutrino event rate is seen at the 2.4σ level. The seasonal modulation of the neutrino flux is also examined, but the change in flux at the SK site is predicted to be negligible, and, as expected, no evidence for a seasonal correlation is seen.