Order of diagenetic events controls evolution of porosity and permeability in carbonates

The properties of carbonate rocks are often the result of multiple, diagenetic events that involve phases of cementation (porosity occlusion) and dissolution (porosity enhancement). This study tests the hypothesis that the order of these events is a major control on final porosity and permeability. A three‐dimensional synthetic model of grainstone is used to quantify trends that show the effect of early cementation, non‐fabric selective dissolution, and then a second‐generation of (post‐dissolution) cement. Models are 3 mm3 with a resolution of 10 μm. Six simple paragenetic sequences are modelled from an identical starting sediment (without accounting for compaction) where the same diagenetic events are placed in different sequences, allowing for quantification of relative changes in the resultant porosity and permeability for each diagenetic event, the trajectory through time, as well as for each final rock. All modelled paragenetic sequences result in reductions in porosity and permeability, but the order of diagenetic events controls the trajectory and final rock properties. Differences in the order of early cement precipitation alone produce variable final values, but all follow the porosity–permeability relationship as expressed by the Kozeny‐Carman equation. However, final values for the sequences which include a phase of dissolution fall on a new curve, which departs from that predicted by the Kozeny‐Carman relationship. This allows an alternative form of porosity–permeability relationship to be proposed: κ = 2280ϕ–30,400, where ϕ is porosity (%) and κ is permeability (mD). Hence while the Kozeny‐Carman relationship predicts porosity–permeability changes that occur with cementation, it is unable to capture accurately changes within the pore network as a result of dissolution. Although the results may be dependent on the properties of the initial carbonate sediment and simplified diagenetic scenarios, it is suggested that this new porosity–permeability relationship may capture some generalized behaviour, which can be tested by modelling further sediment types and diagenetic scenarios.


INTRODUCTION
Carbonate rocks are very heterogeneous, due to both the variability of sediment grain types and the complexity of subsequent diagenetic overprint. From a single sediment, many paragenetic sequences composed of differing successive diagenetic events are possible, leading to a multitude of final rock expressions with differing properties. In carbonates, early cementation tends to occur rapidly and before significant mechanical compaction but may be followed by phases of dissolution and later cementation. While changes in porosity with cementation and dissolution can be readily quantified, the corresponding change in permeability is far more complex and difficult to predict because it is controlled by the precise location of precipitation or dissolution at the pore scale. For example, chemical compaction closes pore throats thus isolating pore bodies and reducing permeability (Budd, 2001). Cementation can also reduce pore throats as well as increasing the surface area and tortuosity of pores, so dramatically reducing permeability. By contrast, dissolution can increase the size of both pores and pore throats leading to enhanced porosity and pore topology (or connectivity) and so increasing permeability.
Empirical porosity-permeability relationships, such as the Kozeny-Carman (KC) equation (Kozeny, 1927;Carman, 1937), are widely used to predict permeability from porosity: where j is the permeability (m 2 ), / is the porosity (fractional), c is the Kozeny constant (m À6 ) and S is the specific surface area (m 2 ). The KC equation assumes, however, uniform changes in porosity with permeability (Carman, 1937). These relationships, therefore, do not allow for either a non-uniform size distribution of pores or pore throats (i.e. heterogeneity) or multiple permeability values for any given porosity value, both of which are the frequent consequence of diagenesis. Pore scale modelling combined with the Lattice Boltzmann method to simulate threedimensional (3D) flow at the pore scale has been widely used to predict carbonate reservoir rock properties such as absolute permeability, capillary pressure and relative permeability. Such models are often based either on statistical reconstruction (see Okabe & Blunt, 2007;Roth et al., 2011;Van der Land et al., 2013) or process, object-based modelling (Biswal et al., 2007;Lander et al., 2008;Mousavi et al., 2012;Harland et al., 2015;Hosa & Wood, 2017). X-ray microtomography (l-CT) can also capture 3D images with a high spatial resolution (ca 1 lm) and reaction rates can be calculated from measured changes in porosity and effluent analysis (Menke et al., 2017). Their general applicability, however, is hindered by the extreme diversity and complexity of natural carbonate rocks. The advantage of synthetic digital rock models is that they are flexible and considerably less expensive than data obtained via laboratory imaging or measurement. They can be used to explore many valuable parameters that relate to the complexity of carbonates, as well as the issues of heterogeneity such as such as Representative Elementary Volumes (REV)the smallest volume over which a measurement can be made that will yield a value representative of the whole. The REVs are in part determined by, and inversely related to, porosity, but also increase with heterogeneity (Menke et al., 2017;Claes et al., 2018). More importantly, due to their flexibility, synthetic digital rock models can establish trends and test hypotheses as to how different depositional and diagenetic factors influence porosity and permeability, as well as to quantify the effect of particular diagenetic events and processes, and indeed the influence of the order of diagenetic events.
In carbonates, early cementation tends to occur rapidly within physically stabilized sediment and before significant mechanical compaction. Early diagenesis can then influence the effect of later diagenetic processes that occur during burial (Van der Land et al., 2013). Here, simple but unique paragenetic sequences (diagenetic pathways) are modelled on a synthetic carbonate sand producing 'digital rocks', that allow exploration of the impact of the order of sequential diagenetic events, including cementation and dissolution, on both intermediate (transient) and final rock properties by altering the paragenetic sequence. This study employs an objectbased model (Hosa & Wood, 2017), an approach which complements studies that simulate diagenesis via reactive transport and reactive pore network models, l-CT imaging and laboratory experiments (De Boever et al., 2012;Menke et al., 2017;Claes et al., 2018). This approach differs from previous object-based models and pore-network simulations (Beckingham, 2017;Hosa & Wood, 2017) in that two types of early marine calcite cements are modelled, as well as dissolution and a later phase of cementation, and then the type and the order of diagenetic events is permutated to produce six different paragenetic sequences. This study does not set out to either accurately model complex diagenetic processes or their products, to validate against real data, or indeed to include all possible diagenetic processes known. Its focus is not on absolute values but rather on relative changes and trends. Simplified diagenetic processes are modelled in order to test the hypothesis that the order of diagenetic events can be a critical control on the final porosity and permeability of a carbonate rock, and to enable quantification of these effects.

METHODS
Carbonate cementation requires fluid oversaturation with respect to carbonate, a source of Ca 2+ and HCO 3 2À such as unstable detrital carbonate, and the presence of suitable nucleation sites.
Parameters that control the geometry of cemented zones include permeability variations that create localized or channelized fluid flow, and the residence time of fluids (see Menke et al., 2016Menke et al., , 2017Claes et al., 2018).
The effect of diagenetic events on a synthetic carbonate sand is modelled using a modified version of Calcite3D ( Fig. 1; Hosa & Wood, 2017). Calcite3D is a process-based and objectbased model with a stochastic initial condition that simulates cementation in natural carbonates, and models two types of early marine calcite cement: syntaxial on monocrystalline grains and isopachous on polycrystalline grains (Fig. 2). Models have a 3 mm 3 volume with a resolution of 10 lm. This sample size and resolution were chosen in order to capture changes D 3D labelled grains E growth of early cement in 2D cross-section  within many pores, and to enable the multiple simulations necessary at reasonable computing cost. The sample size in this study exceeds the proposed volume of REVs shown via numerical simulations to capture heterogeneity due to calcite cement dissolution, from 1.13 mm 3 in low porosity, to 0.2 mm 3 in high porosity, sandstones (Claes et al., 2018). In addition, such a sample size makes the models herein broadly comparable to data derived from l-CT imaging and experiments (Menke et al., 2017).
Diagenetic events are simulated upon uncompacted sediments. Syntaxial cement is modelled as a blocky rhombohedral form with relatively low aspect ratio ( Fig. 2A). Isopachous cement forms fringes around grains (Fig. 2B). Calcite3D takes into account the geometries of the evolving cements, as well as local effects such as grain impingement (Fig. 2C). Isopachous cements are often the earliest cements formed in the marine phreatic realm, but syntaxial cements have been noted to grow synchronously in Bermuda beachrock (Perkins, 1985). Here, Calcite3D is expanded to include two further diagenetic processes: non-fabric-selective dissolution and a phase of later cementation.
The Lattice Boltzmann method is utilized on the entire 3 mm 3 volume to quantify permeability with an overall goal to establish trends in the evolution of porosity and permeability to different diagenetic events. Permeability was quantified for a subset (19 runs for each of the six modelled paragenetic sequences) of fully cemented model outputs on the Archer UK National Supercomputer. Additionally, for several samples, permeability was quantified through all stages in cementation (after each one voxel thick layer of cement was added during one iteration in Calcite3D). Permeability was simulated in three orthogonal directions (x, y and z) through synthetic samples from which the mean was derived. Pore size distribution and connectivity information was extracted using Pore Architecture Tools (PATs; Jiang et al., 2007Jiang et al., , 2012. In the case of monocrystalline grains, the orientation of the crystal axis of the grain is selected at random and the bounding polyhedron is determined ( Fig. 2A). This bounding polyhedron outlines the syntaxial cement growth that will occur in a fast epitaxial stage of cement growth on non-euhedral faces of the grain. The growth of syntaxial cement is modelled until the fast, epitaxial stage of growth is completed.
Syntaxial cement on every monocrystalline grain grows until the bounding polyhedron for that grain is filled completely, or stops earlier if the growth is obstructed by impinging grains. Isopachous cement growth on polycrystalline grains grows up to a thickness controlled by input parameter b.
Later, post-dissolution, cement growth is also modelled by a further isopachous cement growth. While this does not closely resemble the equant crystal form or centripetal growth style of late burial calcite spar, isopachous cement growth is used to approximate the common bladed prismatic form of early burial cement that can grow on earlier cements. Dissolution is a one-step process based on the flow velocity field in the pore space. After the flow field is obtained in the cemented rock via Lattice Boltzmann simulation, this is used to calculate the 'dissolving radius' for each voxel in the pore space, which is a function of the flow velocity at that voxel. The function in the model herein is logarithmic and has the following form: where r d is the radius of dissolution for a given pore voxel, |v| is the magnitude of velocity of flow in that voxel in lattice units as simulated in one direction through the cemented grainstone model output and c is a constant. Integer values of c from 10 to 14 were interrogated, which had a significant impact on dissolution; c = 12 was chosen as this was the value that most consistently removed ca 10% of the solid phase. This functional formulation captures the established concept that if there is higher flow within a channel, more material will be removed (Menke et al., 2016(Menke et al., , 2017. A logarithmic function to dampen the effect of very high velocity flow, which usually occurs in the middle of the channel, with relation to the pore/solid interface. Dissolution is implemented as removal of the solid material by turning all solid voxels inside the radius of dissolution of any pore into pore. The same function is used to calculate the 'dissolving radius' of pores in all simulated paragenetic sequences. Synthetic grainstones with any number of polycrystalline and monocrystalline grains can be generated, but here a single sediment sample is generated with ca 50% monocrystalline grains and used to investigate the evolution of the porosity and permeability for six different hypothetical paragenetic sequences (Fig. 3). Grain size for both grain types follows a normal distribution from 0.05 to 1.3 mm, with a mean of 0.25 mm, and comprises a mixture of spherical to elliptical shapes (Fig. 3A).
All sequences start with a phase of cementation (C), which can involve the alternation of syntaxial (s) and isopachous (i) early marine cement and their growth in parallel (p). Three pathways (CisD, CpD and CsiD) end with a phase of dissolution (D), and three pathways (CisDC, CpDC and CsiDC) follow similar early cementation events, then a phase of dissolution, but end with a further phase of simulated later isopachous cementation (C).
In the parallel mode of cementation one voxel layer of cement is added in each iteration first around all polycrystalline grains (isopachous cement) and then around all monocrystalline grains (syntaxial cement is only added within the limits of the bounding polyhedron). In syntaxialfirst mode, all monocrystalline grains develop syntaxial cement first, and only after the last grain fills its bounding polyhedron completely is isopachous cement allowed to grow on the polycrystalline grains. The number of steps necessary to complete the syntaxial cement growth are unknown a priori, because they depend on the size and orientation of the bounding polyhedra in the synthetic sample, as well as on the neighbourhood of any given grain and whether its growth is obstructed by impinging grains. In isopachous-first mode, all isopachous cement is developed first. The number of iterations for this to occur is fixed by the input parameter b. After the isopachous cementation stage is completed syntaxial cement growth follows on the monocrystalline grains.
In order to keep all modelling parameters the same across all six modelled pathways, two steps of isopachous cement growth are implemented where the volume of the first phase of cementation is reduced in pathways where there is also a second phase of cementation (CisDC, CpDC and CsiDC). Early isopachous cement is modelled as 30 lm thick and later isopachous cement as 20 lm thick on polycrystalline grains that post-dates dissolution. In this way the cumulative thickness of isopachous cement growth (50 lm) is identical for all simulated paragenetic sequences (Fig. 3).
Permeability was quantified via Lattice Boltzmann simulations for the final expressions of synthetic rocks as well as throughout the entire simulated diagenesis for the six modelled pathways. Changing connectivity is explored, as expressed by Euler number, which is the number of nodes (or pores) minus the number of bonds (or pore throats) divided by the volume of the model output. The Euler number is used to derive the pore size dependent connectivity function by removing pores from the network in order of increasing size and recalculating the Euler number at each step (Vogel & Roth, 2001). Removing progressively larger pores will eventually result in a positive Euler number (fewer bonds than nodes), which indicates the point at which the network becomes unconnected. As such, an Euler number of zero is a proxy for the percolation radius of the network.
Outputs are compared to the porosity-permeability (poroperm) relationship as expressed by Global Hydraulic Elements curves (GHE) (Corbett & Potter, 2004), defined by curves based on the KC equation and have the following form: where j is the permeability (mD), / is the porosity (fractional) and FZI is Flow Zone Indicator, often used in reservoir rock typing. The FZI is an absolute measure of pore-throat size in relation to porosity, where a bigger FZI number denotes larger pore throats in relation to pores (m,FZI values 192,96,48,24,12,6 and 3 are used for plots in this paper). In the GHE framework, a regular series of FZI values (an exponential series in permeability space) based on a simple progression is selected to define different fields.

RESULTS
The initial sediment has a porosity of 37.8% and permeability of 57D. All modelled paragenetic sequences result in reductions in porosity (ca 30 to 16%), as well as in permeability (35 to 8D) (Figs 4 and 5). Flow simulations from the three different orthogonal flow directions show very similar poroperm values compared to the mean flow (Fig. 4A).
Early cementation alone occludes pore space and so all three scenarios modelled show a decrease of porosity and permeability (Fig. 4). In pathways where the final step is dissolution (C*D), the isopachous-first cementation mode  (CisD), has a higher volume of isopachous cement and retains the highest porosity and permeability (13.4% and 31D, respectively; Figs 4A and 6A). While the parallel cementation mode (CpD) has less syntaxial cement than the sample with syntaxial-first cementation, it has more than the sample with the isopachous-first cement and shows an intermediate porosity and permeability reduction (17.0% and 39D, respectively; Figs 4B and 6B). The syntaxial-first cementation mode (CsiD), has a higher volume of syntaxial cement than isopachous cement and shows the most marked decrease in porosity and permeability, to 20.4% and 47D, respectively (Figs 4C and 6C).
Pathways where the final step is dissolution (C*D; Fig. 5A to C) show lower final porosities (17.5 to 24.4%) and permeabilities (10 to 26D), than pathways where the final step is cementation (C*DC; Fig. 5D to F), which show porosities of 24.8 to 29.7%, and permeabilities of 25 to 40D.
Notably, a single curve is defined by the final poroperm values of all simulations (j = 2276/-30429, simplified to three significant numbers to j = 2280/-30400, where / is porosity and j is permeability) and can be fitted with a linear function with a high coefficient of determination (R 2 = 0.997) (Fig. 5A). The FZI number for the six model end-points is also a linear function of both porosity (R 2 = 0.998; Fig. 5B) and permeability (R 2 = 0.991; Fig. 5C). Early cementation trends alone in all pathways follow the porosity-permeability relationship as expressed by GHE curves and the KC equation (Fig. 4). The phase of dissolution, however, notably moves the evolving rock to another FZI level, changing the FZI number by between 3 to 18 units.
Pore size distributions in all models are simple and effectively unimodal (Fig. 7), but pore size distribution changes dramatically through the six pathways ( Fig. 7A to C). The three pathways where the final step is cementation (C*DC) have a higher proportion of larger pores remaining compared to C*D pathways, even though the absolute number of large pores is lower (Fig. 7A to C). Likewise, connectivity, as expressed by Euler number, in the C*DC pathways results in a higher connectivity between the larger pores (>100 lm diameter) than C*D pathways (Fig. 7D to F).

Limitations of the model
One initial sediment grain-size distribution is used, and no other grain shapes apart from near spheroids are employed. The effectively unimodal pore size distributions of all models (Fig. 7) correspond to only some natural carbonates, such as the Portland limestone commonly used in experiments (Menke et al., 2017). The simulated porosities and permeabilities are also far higher than natural carbonate rocks. Initial porosity and permeability values of the sediment are 37.8% and 57MD. While these are comparable to modern grainy, unconsolidated, carbonate sand, they do not account for any post-depositional mechanical or chemical compaction in the models. All diagenetic events in this study are therefore simulated upon uncompacted sediments, and so do not correspond to any types of 'late' diagenesis which can occur in either the burial realm or the near surface as a result of uplift (i.e. telogenesis).
The models of diagenetic processes are also simplified, for example, the last phase of cementation is grown isopachously to simulate bladed prismatic early burial cement. This does not approximate to late burial equant cements that fill pores from pore edges inwards by the growth of a limited number of large crystals, often leading to the partial or complete occlusion of macropores. Burial cementation is also not dependent upon a prior isopachous cement, as has been modelled here.
Likewise the dissolution of grains to form moulds is not considered, nor any other constants of dissolution interrogated or variable quantiles of material removed other than 10%. For example, a different poroperm relationship was noted in the experimental dissolution and cementation using l-CT imaging and numerical modelling of dolomite/evaporitic rocks (De Boever et al., 2012). Fractures, which are common in carbonate rocks and often impart substantial permeability because they provide the largest, best-connected pores, are not considered in this study.
The models also have a relatively coarse resolution (10 lm), representing a trade-off between achieving a sufficient volume of sample interrogated, with computational cost. At such a resolution, no features such as asperities due to finescale cementation or dissolution, or the contribution of micropores to flow, can be captured.
However, the intention here is to test the simple hypothesis: does the order of diagenetic events matter? Given that all final simulations fall on the same curve, it is suggested that a broadly general set of behaviours may have been identified. Yet other pore topologies, and perhaps different results, might arise if other scenarios are tested. The validity and predictability of all models depends on the representativeness of the pore networks and on the equations and parameters used to model the diagenetic events.

Quantifying the effects of diagenesis
Simulation results show that in all paragenetic sequences, early cementation alone decreases porosity, but permeability decreases more dramatically and variably. All cementation simulations follow that predicted by the standard Kozeny-Carman porosity-permeability relationship.
While early cementation trends in all pathways closely follow the porosity-permeability relationship as expressed in the KC equation, dissolution, by contrast, moves the evolving rock to another KC curve (Fig. 6). Dissolution is therefore controlled by the pre-existing pore topology, which is then markedly altered by dissolution to create a new pore network. Indeed more dissolution occurs, in the C*DC samples than in the C*D, and this extra dissolution is what determines the ultimate outcome, that C*D has higher final poroperm.
That dissolution does not follow the predicted porosity-permeability relationship of KC equations has been found in pore-network model simulations (Beckingham, 2017). Kozeny-Carman relationships appear to be unable to accurately capture the changes occurring in the pore network given this type of diagenetic reaction. A similar conclusion was also reached using dynamic reservoir-condition microtomography of reactive transport in complex carbonates  (Menke et al., 2017), which found that in addition to brine pH, it is pore structure that determines the style of dissolution. Where porosity is homogenous, dissolution is uniform, but where the porosity increase is concentrated in preferential flow paths, i.e. is heterogeneously distributed, dissolution occurs as channelling and, in increasingly complex carbonates, channelized flow forms faster.
All simulations fall on a single curve, which again does not follow a KC relationship. This is probably a consequence of the variable way in which dissolution affects the rock, controlled by pre-dissolution properties, which in turn is dictated by early cementation, which differs for the six paragenetic sequences. The form of this relationship is likely to be dependent on the initial sediment grain-size distribution, as increasing grain size will increase permeability but not affect porosity.
Dissolution has the notable effect of changing the FZI (DFZI). Pre-dissolution porosity (Fig. 5D), permeability (Fig. 5E) and FZI (Fig. 5F) control the change in FZI during dissolution in the models. More dissolution results in an increased shift in DFZI pre-dissolution and post-dissolution (Fig. 5D), and even small differences in FZI predissolution can lead to significant differences in FZI post-dissolution (range 2, Fig. 5F), as well as porosity (Fig. 5D) and permeability (Fig. 5E). This observation has implications for the determination of petrophysical reservoir rock types.
Results from this study demonstrate that the order of diagenetic events has a marked impact on final carbonate porosity and permeability. Counterintuitively, the three pathways that end with a phase of later cementation (C*DC) result in higher porosities and permeabilities, an overall greater pore size, and higher connectivity between the larger pores than the three pathways that end with dissolution (C*D). It is likely that cementation fully occludes smaller pores first such that the proportion of larger pores remaining increases even though the absolute number of large pores is lower. When the Euler number connectivity curve passes though zero, this indicates the pore size at which the pore network has become disconnected. The three C*D pathways disconnect at ca 50 lm, but the C*DC pathways disconnect at ca 60 lm (Fig. 7F). This suggests that medium-sized pores (ca 50 lm) are better connected in the C*DC pathways. This is likely due to the smaller volume of early cementation in C*DC, such that the ca 50 lm pores were not cut-off from one another before dissolution. In the dissolution model, once the pores are disconnected, there is minimal chance of re-connection during dissolution, since the dissolution model only affects the solid/pore interface where there is fluid flow.
In the Cis*, Cp* and Csi* pathways, 50 lm pores are still connected for Cis*, but disconnected for Csi* (Fig. 7F). As before, in the two Csi* pathways the volume of pre-dissolution cement is the greatest (Fig. 4C and F), so the pores are cut-off from one another to a greater extent than in Cis* ( Fig. 4A and D), thus resulting in greater disconnection of pores after dissolution.
Although the models do not quantify changes in pore-throat size, it can be hypothesized that a reduction in pore-throat size and hence loss of flow-conducting pores is responsible for reducing permeability. So, while the KC equation predicts a reduction in permeability, this simple relationship between porosity and permeability does not account for the effects of changes in pore-throat sizes (Crandell et al., 2012). An alternative form of porosity-permeability curve that captures the end-points of the modelled processes simulated via six paragenetic sequences is therefore proposed: j = 2280/-30400, where / is porosity in % and j is permeability in mD (Fig. 5A).
An investigation of the relationship between porosity and permeability and the FZI values also reveals that these can be described as linear, i.e. with high coefficients of determination R 2 (Fig. 5D to F). The change in FZI during dissolution, which can be viewed as a proxy for pore topology, is controlled by porosity and permeability. The lower the pre-dissolution poroperm, the bigger the change in FZI during dissolution ( Fig. 5D to F). The variation in FZI values during cementation does not exceed 6.2 for any of the paragenetic sequences (reaching values 6.2, 4.1, 3.1, 3.8, 4.1 and 2.3 for CisD, CpD, CsiD, CisDC, CpDC and CsiDC, respectively; Fig. 5). However, the change in FZI during dissolution (ΔFZI DISSOLUTION ) can be significantly higher, especially for the pathways with no later cementation (C*D), where this can reach up to 17.4 (reaching values 10.5, 14.3, 17.4, 3.0, 3.5 and 6.5 for CisD, CpD, CsiD, CisDC, CpDC and CsiDC, respectively; Fig. 5).
The final FZI will also be biggest for the smallest pre-dissolution poroperm. Since the pre-dissolution poroperm is controlled by the mode and volume of early cementation, it can be concluded that early cementation also plays a role in controlling the final pore topology.
Due to the convex shape of the early cementation curves (Fig. 4), ultimately a vertical slope would be achieved with continued cementation, so this suggests that there is a critical porosity when precipitation begins to impact the smallest pores, and permeability decreases rapidly. By extrapolating the primary cementation curves with an exponential function (j = ae Àb/ + c, where j is the permeability [mD], / is the porosity [%]), the following values of porosity are obtained where permeability is zero: CisDC: 13.3%; CsiDC: 13.4%; CpDC: 13.1%; CisD: 11.4%; CsiD: 11.2%; CpD: 9.7%. While variation in the size of the smallest pores has little impact on porosity in these models, this has a significant impact on permeability. Predictions from pore-network modelling suggest that this is controlled by the smallest conducting pore-throat size as they are often highly abundant (Beckingham et al., 2013). Small changes in pore-throat size due to diagenesis can then have a marked effect on permeability. This is supported by the experimental dissolution of a dolomite/evaporitic rock, where initial dissolution led to an overall increase in pore size due to pore rather than pore-throat dissolution, but further dissolution resulted in the loss of the largest pores, but a marked increase in permeability due to enlargement of pore throats (De Boever et al., 2012).

CONCLUSIONS
Models of different diagenetic scenarios demonstrate that the order of events has a notable control on both intermediate and final carbonate rock properties. These also demonstrate that current porosity-permeability relationships such as the Kozeny-Carman (KC) equation are unable to reflect the pore scale phenomena that control these relationships. This is because diagenesis creates size-dependent reactions which, critically, infer that changes in the smallest conducting pore throats can have a marked effect on permeability. Differences between the styles of early cementation also exert a major control on the evolution of poroperm, because this determines the pore topology prior to dissolution, which in turn controls how dissolution will affect the evolution of rock properties. Moreover, while cementation alone tends to follow the KC curve, dissolution has a notable effect by changing this trajectory. On the basis of these models an alternative form of porosity-permeability relationship curve is proposed: j = 2280/-30400, where / is porosity in % and j is permeability in mD.
Although results in this study may be dependent on the properties of the initial carbonate sediment and simplified diagenetic scenarios, it is suggested that this new porosity-permeability relationship may capture some generalized behaviour, which can be tested by modelling additional and more complex diagenetic scenarios.