LIGO and Virgos third observing run (O3) revealed the first neutron star-black hole (NSBH) merger candidates in gravitational waves. These events are predicted to synthesize r-process elements creating optical/near-IR kilonova (KN) emission. The joint gravitational-wave (GW) and electromagnetic detection of an NSBH merger could be used to constrain the equation of state of dense nuclear matter, and independently measure the local expansion rate of the universe. Here, we present the optical follow-up and analysis of two of the only three high-significance NSBH merger candidates detected to date, S200105ae and S200115j, with the Zwicky Transient Facility (ZTF). ZTF observed $sim$,48% of S200105ae and $sim$,22% of S200115js localization probabilities, with observations sensitive to KNe brighter than $-$17.5,mag fading at 0.5,mag/day in g- and r-bands; extensive searches and systematic follow-up of candidates did not yield a viable counterpart. We present state-of-the-art KN models tailored to NSBH systems that place constraints on the ejecta properties of these NSBH mergers. We show that with depths of $rm m_{rm AB}approx 22$ mag, attainable in meter-class, wide field-of-view survey instruments, strong constraints on ejecta mass are possible, with the potential to rule out low mass ratios, high BH spins, and large neutron star radii.