The results of a search for Higgs boson pair () production in the , , and decay modes are presented. The search uses proton-proton collision data recorded by the CMS experiment at a center-of-mass energy of 13~\mathrm{TeV}, corresponding to an integrated luminosity of 138~\mathrm{fb}^{-1}. The analyzed events contain two, three, or four reconstructed leptons, including electrons, muons, and hadronically decaying tau leptons. No evidence for a signal is found in the data, and upper limits are set on the cross section for nonresonant production, as well as resonant production in which a new heavy particle decays to a pair of Higgs bosons. For nonresonant production, the observed (expected) upper limit on the cross section at confidence level (CL) is () times the standard model (SM) prediction. The ratio of the trilinear Higgs boson self-coupling to its value in the SM is constrained to be between and at CL, and limits are set on a variety of new-physics models using an effective field theory approach. The observed (expected) limits on the cross section for resonant production range from to ( to )~\mathrm{pb} at CL for new heavy-particle masses from to 1000~\mathrm{GeV}.