A search for supersymmetry (SUSY) is presented where at least one Higgs boson is produced and decays to two photons in the decay chains of pair-produced SUSY particles. Two complementary analysis strategies are pursued: one focused on strong SUSY production and the other focused on electroweak SUSY production. The presence of charged leptons, additional Higgs boson candidates, and various kinematic variables are used to categorize events into search regions that are sensitive to different SUSY scenarios. The results are based on data from proton-proton collisions at the Large Hadron Collider at a center-of-mass energy of 13 TeV collected by the CMS experiment, corresponding to an integrated luminosity of 77.5 fb$^{-1}$. No statistically significant excess of events is observed relative to the standard model expectations. We exclude bottom squark pair production for bottom squark masses below 530 GeV and a lightest SUSY particle mass of 1 GeV; wino-like chargino-neutralino production for chargino and neutralino masses below 235 GeV with a gravitino mass of 1 GeV; and higgsino-like chargino-neutralino production in the case where the neutralino decays exclusively to a Higgs boson and a gravitino for neutralino masses below 290 GeV.