  By employing Monte Carlo random sampling, traditional binary population
synthesis (BPS) offers a substantial improvement in efficiency over brute
force, grid-based studies. Even so, BPS models typically require a large number
of simulation realizations, a computationally expensive endeavor, to generate
statistically robust results. Recent advances in statistical methods have led
us to revisit the traditional approach to BPS. In this work we describe our
publicly available code dart_board which combines rapid binary evolution codes,
typically used in traditional BPS, with modern Markov chain Monte Carlo
methods. dart_board takes a novel approach that treats the initial binary
parameters and the supernova kick vector as model parameters. This formulation
has several advantages, including the ability to model either populations of
systems or individual binaries, the natural inclusion of observational
uncertainties, and the flexible addition of new constraints which are
problematic to include using traditional BPS. After testing our code with mock
systems, we demonstrate the flexibility of dart_board by applying it to three
examples: (i) a generic population of high mass X-ray binaries (HMXBs), (ii)
the population of HMXBs in the Large Magellanic Cloud (LMC) in which the
spatially resolved star formation history is used as a prior, and (iii) one
particular HMXB in the LMC, Swift J0513.4-6547, in which we include
observations of the system's component masses and orbital period. Although this
work focuses on HMXBs, dart_board can be applied to a variety of stellar
binaries including the recent detections by gravitational wave observatories of
merging compact object binaries.
