Simulating low-likelihood high-impact extreme weather events in a warming world is a significant and challenging task for current ensemble forecasting systems. While these systems presently use up to 100 members, larger ensembles could enrich the sampling of internal variability. They may capture the long tails associated with climate hazards better than traditional ensemble sizes. Due to computational constraints, it is simply impossible to generate huge ensembles, of say 10,000 members, using traditional numerical simulations of climate models at high resolution. We replace traditional numerical simulations with machine learning (ML) models to generate hindcasts of huge ensembles. We construct an ensemble weather forecasting system based on Spherical Fourier Neural Operators, and we discuss important design decisions for constructing such an ensemble. The ensemble represents model uncertainty through perturbed-parameter techniques, and it represents initial condition uncertainty through bred vectors, which sample the fastest growing modes of the forecast. Using the IFS operational weather forecasting system as a baseline, we present a rigorous diagnostics pipeline. We evaluate the ML ensemble's overall performance as a probabilistic forecasting system and show that its performance is comparable to (though approximately 12-18 hours behind) that of IFS. Based on the spectra of the individual ensemble members and the ensemble mean, the ML ensemble trajectories have realistic error growth. We also evaluate the extreme weather forecasts in particular, by assessing the ensemble's reliability and discrimination. These diagnostics test the physical fidelity of the ML emulator. They ensure that the ensemble can reliably simulate the overall time evolution of the atmosphere, including low likelihood high-impact extremes. We generate a huge ensemble initialized each day in summer 2023, and we characterize the statistics of extremes.