Many astrophysical systems, such as stellar binaries, require detailed numerical simulations to model. To explore the effects of uncertain parameters on final outcomes, we often use regularly spaced grids of simulations where the density and spacings are arbitrarily chosen using some combination of intuition, visual inspection, and computational cost per simulation. This grid-based approach becomes more inefficient and intractable for high dimensional studies. Moreover, this approach falls short for studies relying on the accuracy of interpolation between simulation results. Binary stellar evolution suffers from both of these problems, existing in a large parameter space and having highly non-linear evolution. We present a machine-learning approach to performing parameter studies using active learning and find improvements in classification and regression when using our algorithm compared to a regularly spaced grid. Our algorithm uses parallel tempered MCMC and combines classification and regression metrics simultaneously for proposing new simulations. We apply it to optimize large grids of binary evolution simulations using the stellar evolution code MESA. Our algorithm will be a key component in POSYDON, the next-generation binary population synthesis code we are actively developing.