The production of Z boson pairs in proton-proton ( p p ) collisions, p p → ( Z / γ ∗ ) ( Z / γ ∗ ) → 2 ℓ 2 ℓ ' , where ℓ , ℓ ' = e or μ , is studied at a center-of-mass energy of 13 TeV with the CMS detector at the CERN LHC. The data sample corresponds to an integrated luminosity of 137 fb - 1 , collected during 2016-2018. The Z Z production cross section, σ tot ( p p → Z Z ) = 17.4 ± 0.3 (stat) ± 0.5 (syst) ± 0.4 (theo) ± 0.3 (lumi) pb , measured for events with two pairs of opposite-sign, same-flavor leptons produced in the mass region 60 < m ℓ + ℓ - < 120 GeV is consistent with standard model predictions. Differential cross sections are also measured and agree with theoretical predictions. The invariant mass distribution of the four-lepton system is used to set limits on anomalous Z Z Z and Z Z γ couplings.