Solutions to the $mu$ problem in supersymmetry based on the Kim-Nilles mechanism naturally feature a Dine-Fischler-Srednicki-Zhitnitsky (DFSZ) axion with decay constant of order the geometric mean of the Planck and TeV scales, consistent with astrophysical limits. We investigate minimal models of this type with two gauge-singlet fields that break a Peccei-Quinn symmetry, and extensions with extra vectorlike quark and lepton supermultiplets consistent with gauge coupling unification. We show that there are many anomaly-free discrete symmetries, depending on the vectorlike matter content, that protect the Peccei-Quinn symmetry to sufficiently high order to solve the strong CP problem. We study the axion couplings in this class of models. Models of this type that are automatically free of the domain wall problem require at least one pair of strongly interacting vectorlike multiplets with mass at the intermediate scale, and predict axion couplings that are greatly enhanced compared to the minimal supersymmetric DFSZ models, putting them within reach of proposed axion searches.