Partially dolomitized carbonate rocks of the Middle East and North America host large hydrocarbon reserves. The origin of some of these dolomites has been attributed to a hydrothermal mechanism. The Benicàssim area (Maestrat Basin, eastern Spain) constitutes an excellent field analogue for fault-controlled stratabound hydrothermal dolomitization: dolostone geobodies are well exposed and extend over several kilometres away from seismic-scale faults. This work investigates the main controls on the formation of stratabound versus massive dolomitization in carbonate sequences by means of two-dimensional (2D) reactive transport models applied to the Benicàssim case study. Simulation results suggest that the dolomitization capacity of Mg-rich fluids reaches a maximum at temperatures around 100 °C and a minimum at 25 °C (studied temperature range: 25–150 °C). It takes of the order of hundreds of thousands to millions of years to completely dolomitize kilometre-long limestone sections, with solutions flowing laterally through strata at velocities of metres per year (m/a). Permeability differences of two orders of magnitude between layers are required to form stratabound dolomitization. The kilometre-long stratabound dolostone geobodies of Benicàssim must have formed under a regime of lateral flux greater than metres per year over about a million years. As long-term dolomitization tends to produce massive dolostone bodies not seen at Benicàssim, the dolomitizing process there must have been limited by the availability of fluid volume or the flow-driving mechanism. Reactive transport simulations have proven a useful tool to quantify aspects of the Benicàssim genetic model of hydrothermal dolomitization.