The cosmic merger rate density of black hole binaries (BHBs) can give us an essential clue to constraining the formation channels of BHBs, in light of current and forthcoming gravitational wave detections. Following a Monte Carlo approach, we couple new population-synthesis models of BHBswith the Illustris cosmological simulation, to study the cosmic history of BHB mergers. We explore six population-synthesis models, varying the prescriptions for supernovae, common envelope and natal kicks. In most considered models, the cosmic BHB merger rate follows the same trend as the cosmic star formation rate. The normalization of the cosmic BHB merger rate strongly depends on the treatment of common envelope and on the distribution of natal kicks. We find that most BHBs merging within LIGO's instrumental horizon come from relatively metal-poor progenitors (<0.2 Z⊙). The total masses of merging BHBs span a large range of values, from ~6 to ~82M⊙. In our fiducial model, merging BHBs consistent with GW150914, GW151226 and GW170104 represent ~6, 3 and 12 per cent of all BHBs merging within the LIGO horizon, respectively. The heavy systems, like GW150914, come from metalpoor progenitors (<0.15 Z⊙). Most GW150914-like systems merging in the local Universe appear to have formed at high redshift, with a long delay time. In contrast, GW151226-like systems form and merge all the way through the cosmic history, from progenitors with a broad range of metallicities. Future detections will be crucial to put constraints on common envelope, on natal kicks, and on the BHB mass function.
