The images in this dataset show the mixing of two liquid solutions in a random bead pack as a function of time and in three-dimensions. The working fluids used in this study are solutions of methanol and ethylene-glycol (MEG, fluid 1) and brine (fluid 2). In particular, three mixtures of ethylene-glycol and methanol were prepared that differ in wt% ethylene-glycol, namely 55 wt% (MEG55), 57 wt% (MEG57) and 59 wt% (MEG59). Measurements are conducted using in the regime of Rayleigh numbers, Ra = 2000-5000. X-ray Computed Tomography is applied to image the spatial and temporal evolution of the solute plume non -invasively. The tomograms are used to compute macroscopic quantities including the rate of dissolution and horizontally averaged concentration profiles, and enable the visualisation of the ow patterns that arise upon mixing at a spatial resolution of about (2x2x2) mm3. We observe that the mixing process evolves systematically through three stages, starting from pure diffusion, followed by convection-dominated and shutdown. A modified diffusion equation is applied to model the convective process with an onset time of convection that compares favourably with literature data and an effective diffusion coefficient that is almost two orders of magnitude larger than the molecular diffusivity of the solute. The comparison of the experimental observations of convective mixing against their numerical counterparts of the purely diffusive scenario enables the estimation of a non-dimensional convective mass flux in terms of the Sherwood number, Sh = 0.025Ra. We observe that the latter scales linearly with Ra, in agreement with observations from both experimental and numerical studies on thermal convection over the same Ra regime.