  Prior to the detection of black holes (BHs) via the gravitational waves (GWs)
they generate at merger, the presence of BHs was inferred in X-ray binaries,
mostly via dynamical measurements, with masses in the range between $\sim
5-20~M_\odot$. The LIGO discovery of the first BHs via GWs was surprising in
that the two BHs that merged had masses of $35.6^{+4.8}_{-3.0}$ and
$30.6^{+3.0}_{-4.4}\,M_\odot$, which are both above the range inferred from
X-ray binaries. With 20 BH detections from the O1/O2 runs, the distribution of
masses remains generally higher than the X-ray inferred one, while the
effective spins are generally lower, suggesting that, at least in part, the
GW-detected population might be of dynamical origin rather than produced by the
common evolution of field binaries. Here we perform high-resolution N-body
simulations of a cluster of isolated BHs with a range of initial mass spectra
and upper mass cut-offs, and study the resulting binary mass spectrum resulting
from the dynamical interactions. Our clusters have properties similar to those
of the massive remnants in an OB association $\sim 10 \, \mathrm{Myr}$ after
formation. We perform a likelihood analysis for each of our dynamically-formed
binary population against the data from the O1 and O2 LIGO/Virgo runs. We find
that an initial mass spectrum $M_{\rm BH}\propto M^{-2.35}$ with an upper mass
cutoff $M_{\rm max}\sim 50M_\odot$ is favored by the data, together with a
slight preference for a merger rate that increases with redshift.
