A measurement is presented of the production of $Z$ bosons that decay into two electrons or muons in association with jets, in proton-proton collisions at a center-of-mass energy of 13 TeV. The data were recorded by the CMS Collaboration at the LHC with an integrated luminosity of $35.9\text{ }\text{ }{\mathrm{fb}}^{\ensuremath{-}1}$. The differential cross sections are measured as a function of the transverse momentum (${p}_{\mathrm{T}}$) of the $Z$ boson and the transverse momentum and rapidities of the five jets with largest ${p}_{\mathrm{T}}$. The jet multiplicity distribution is measured for up to eight jets. The hadronic activity in the events is estimated using the scalar sum of the ${p}_{\mathrm{T}}$ of all the jets. All measurements are unfolded to the stable-particle level and compared with predictions from various Monte Carlo event generators, as well as with expectations at leading and next-to-leading orders in perturbative quantum chromodynamics.