Multiple populations are ubiquitous in the old massive globular clusters (GCs) of the Milky Way. It is still unclear how they arose during the formation of a GC. The topic of iron and metallicity variations has recently attracted attention with the measurement of iron variations among the primordial population (P1) stars of Galactic GCs. We use the spectra of more than 8000 RGB stars in 21 Galactic GCs observed with MUSE to derive individual stellar metallicities [M/H]. For each cluster, we use the HST photometric catalogs to separate the stars into two main populations (P1 and P2). We measure the metallicity spread within the primordial population of each cluster by combining our metallicity measurements with the stars {Delta}_F275W,F814W_ pseudo-color. We also derive metallicity dispersions ({sigma}_[M/H]_) for the P1 and P2 stars of each GC. In all but three GCs, we measure a significant correlation between the metallicity and the {Delta}_F275W,F814W_ pseudo-color of the P1 stars such that stars with larger {Delta_F275W,F814W_ have higher metallicities. We measure metallicity spreads that range from 0.03 to 0.24dex and correlate with the GC masses. As for the intrinsic metallicity dispersions, when combining the P1 and P2 stars, we measure values ranging from 0.02 dex to 0.08dex that correlate very well with the GC masses. We compared the metallicity dispersion among the P1 and P2 stars and found that the P2 stars have metallicity dispersions that are smaller or equal to that of the P1 stars. We find that both the metallicity spreads of the P1 stars (from the {Delta_F275W,F814W_ spread in the chromosome maps) and the metallicity dispersions ({sigma_[M/H]_) correlate with the GC masses, as predicted by some theoretical self-enrichment models presented in the literature.