Gluon splitting to b-quark pairs is shown to be an important process which is currently not well modelled by Monte Carlo event generators. In order to understand the splitting process in Monte Carlo generation, Pythia 8 and Sherpa are compared in terms of their production methods of b-quark pairs, and the parton shower approximation of gluon splitting is compared to its matrix element calculation using Sherpa. The parton shower in Sherpa shows reasonable agreement with the matrix element, but large differences are seen between the predictions of Pythia 8 and Sherpa, particularly in their b-quark pair production through the parton shower. In order to compare the generators to data, a non-prompt di-J/psi sample is prepared via a four-dimensional simultaneous binned-fit to data collected by the ATLAS detector in the Large Hadron Collider at CERN. The Monte Carlo result distributions are then smeared based on B->J/psi correlation functions derived from the analysis of a Pythia 8B J/psi sample. The generator and collider samples are then compared in terms of the distribution of the angular separations between J/psi pairs. It is found that of the two generators: Pythia 8 shows the best agreement with the data, and Sherpa primarily disagrees in the region of low angular-separation, where parton shower production is seen to dominate