We propose and test a hierarchical high-throughput screening (HHTS) approach to catalyst design for complex catalytic reaction systems that is based on quantum mechanics (QM) derived full reaction networks with QM rate constants but simplified to examine only the reaction steps likely to be rate determining. We illustrate this approach by applying it to determine the optimum dopants (our of 35 candidates) to improve the turnover frequency (TOF) for the Fe-based Haber-Bosch ammonia synthesis process. We start from the QM-based free-energy reaction network for this reaction over Fe(111), which contains the 26 most important surface configurations and 17 transition states at operating conditions of temperature and pressure, from which we select the key reaction steps that might become rate determining for the alloy. These are arranged hierarchically by decreasing free-energy reaction barriers. We then extract from the full reaction network, a reduced set of reaction rates required to quickly predict the effect of the catalyst changes on each barrier. This allows us to test new candidates with only 1% of the effort for a full calculation. Thus, we were able to quickly screen 34 candidate dopants to select a small subset (Rh, Pt, Pd, Cu) that satisfy all criteria, including stability. Then from these four candidates expected to increase the TOF for NH3 production, we selected the best candidate (Rh) for a more complete free-energy and kinetic analysis (10 times the effort for HHTS but still 10% of the effort for a complete analysis of the full reaction network). We predict that Rh doping of Fe will increase the TOF for NH3 synthesis by a factor of similar to 3.3 times compared to Fe(111), in excellent agreement with our HHTS predictions, validating this approach.