To study metamorphic carbonation at greenschist facies conditions in the SW Scottish Highlands, a novel inverse modeling framework, which combines solutions of the transport equation with a global optimization method of differential evolution, was developed. Using this framework, we calculated simultaneously time-integrated and time-averaged metamorphic fluid fluxes of 83.4 +/- 35.4 m(3) m(-2) and 10(-10.1) (+/-) (0.5) m(3) m(-2) s(-1), respectively, a time-averaged reaction rate constant of 10(-12.7--10.2) m s(-1) and comparable timescales for fluid flow and fluid-driven reaction of 10(4.3 +/- 0.5) yr and 10(2.7-5.2) yr, respectively. These parameters were calculated using an empirical estimate of the coefficient of molecular diffusion and a calculated value for metamorphic porosity. Our estimates are (1) consistent with single pass flow of fluid released by metamorphic devolatilization, (2) within the range where heat is transported by conduction and matter is transported by advection, (3) in agreement with an emerging consensus that metamorphic events are relatively short-lived, and (4) supportive of applying laboratory-based estimates of kinetic parameters to metamorphic systems. Based on a sensitivity analysis, we show that (1) selecting the diffusion coefficient (rather than fluid velocity, reaction rate or flow duration) as an input parameter yields more robust estimates of metamorphic fluid flow parameters, and (2) ignoring reaction-dependent porosity and reaction rates can result in an order-of-magnitude uncertainty in best-fit flow parameters, evaluated from concentration profiles. Finally, similarity between our calculated time-averaged metamorphic fluid fluxes which were obtained numerically and those which were obtained analytically confirms the validity of using the 'quasi-stationary state' assumption to quantify metamorphic fluid flow parameters.