Impact of textural patterns on modeled rock weathering rates and size distribution of weathered grains

Rock texture has a critical influence on the way rocks weather. The most important textural factors affecting weathering are grain size and the presence of cracks and stylolites. These discontinuities operate as planes of mechanical weakness at which chemical weathering is enhanced. However, it is unclear how different rock textures impact weathering rates and the size of weathered grains. Here, we use a numerical model to simulate weathering of rocks possessing grain boundaries, cracks, and stylolites. We ran simulations with either synthetic or natural patterns of discontinuities. We found that for all patterns, weathering rates increase with discontinuity density. When the density was <~25%, the weathering rate of synthetic patterns followed the order: grid > honeycomb > Voronoi > brick wall. For higher values, all weathering rates were similar. We also found that weathering rates decreased as the tortuosity of the pattern increased. Moreover, we show that textural patterns strongly impact the size distributions of detached grains. Rocks with an initial monomodal grain size distribution produce weathered fragments that are normally distributed. In contrast, rocks with an initial log‐normal size distribution produce weathered grains that are log‐normally distributed. For the natural patterns, weathering produced lower modality distributions.


| INTRODUCTION
Both natural and anthropogenic processes are affected by the rate at which rocks weather. Weathering rates impact the development of landscapes, the formation of soils, the fluid flow in aquifers and petroleum reservoirs, and the durability of buildings and monuments (Brantley, 2008;Dixon et al., 2012;Wilson, 2004). In addition, weathering rate plays a significant role in the global carbon cycle (Li & Elderfield, 2013;Torres et al., 2016), regulating atmospheric CO 2 on geological time scales (Them et al., 2017). Artificially accelerated weathering has even been suggested as a way of mitigating presentday anthropogenic carbon emissions (Beerling et al., 2020;Strefler et al., 2018;Torres et al., 2016;Xu & Liu, 2010).
At the microscopic scale, weathering rates are affected by discontinuities that include crystalline defects, crystal edges and corners, and grain boundaries (Holdren & Speyer, 1987;Trindade Pedrosa et al., 2019). For example, the rate of dissolution along the edges and corners of a calcite spar was measured to be 1.7-3.6 faster than on the mineral face (Noiriel et al., 2019). In polycrystalline rocks, grain boundaries were found to be an order of magnitude more reactive than the bulk mineral (Bray et al., 2015;Emmanuel, 2014;Jonas et al., 2014). In studies focused on rock weathering at the submicron scale, enhanced dissolution at grain boundaries was shown to cause the mechanical detachment of particles into the fluid phase (Emmanuel & Levenson, 2014;Fischer & Luttge, 2017;Krklec et al., 2016;Silveira & Aarão Reis, 2013). Such chemo-mechanical rock weathering was observed in micritic limestone (Levenson & Emmanuel, 2016), however particle detachment can occur in various types of rocks with larger grain sizes and different mineral compositions (Israeli & Emmanuel, 2018;Krklec et al., 2013;Levenson & Emmanuel, 2016).
At macroscopic scales, rock weathering is accelerated by other types of discontinuities, such as cracks, joints, fractures, and stylolites (Heap et al., 2018;Singhal & Gupta, 2010), which operate as planes of mechanical weakness and enhanced chemical weathering (Eppes & Keanini, 2017;Lei et al., 2017;Pacheco & Alencoão, 2006). For example, Røyne et al. (2008) showed that outcrop weathering is controlled by continual fracturing and production of surface area, which allows fluids to penetrate deeper into the rock and accelerate weathering rates.
Furthermore, similar patterns can have different levels of connectivity depending on the spacing, orientation, length, and density of the discontinuities. The convolution of these factors can be represented by tortuosity, which is a measure of the geometric complexity of the pathways by which reactive fluids penetrate the rock. High tortuosity is expected to lead to reduced weathering rates, while low tortuosity could intensify weathering.
Here, we develop a cellular automaton model that simulates coupled chemo-mechanical weathering processes of rocks with different kinds of discontinuities and textural patterns. Specifically, we analyze the impact of the density and tortuosity of the discontinuities on the weathering rate. In addition, we examine how these parameters impact the size distribution of weathered rock fragments. We also discuss the implications for both surface and subsurface processes including soil and regolith production.

| Model structure
To simulate the effects of chemo-mechanical weathering on rocks with different textures and grain size distributions, we used a model based on that described by Israeli & Emmanuel (2018). A 2D crosssection of the rock was represented using a domain with 560 × 420 elements. Each element represented either a solid mineral, a discontinuity, or a fluid phase and is assigned a characteristic value.
These values reflect the element's relative kinetic rate of dissolution.
In our configuration, fluid phase elements were given a characteristic value of zero, discontinuities were given the maximum characteristic value (1), and solid minerals were given a value of 0.1, an order of magnitude smaller than discontinuities, which is consistent with literature values (Bray et al., 2015;Emmanuel, 2014;Jonas et al., 2014).
In the simulations, chemical weathering only occurs in elements neighboring the fluid phase. In every time step, the probability that an element will dissolve depends both on the characteristic value of the element and the number of neighboring fluid elements. The dissolved elements are then reassigned as a fluid phase, and the domain is scanned for interconnected elements that are fully surrounded by fluid. These surrounded elements are considered to be detached physically and their elements are also reassigned to the fluid phase ( Figure 1; see also Supplementary Video S1). Israeli & Emmanuel (2018) used the model with different ratios of reactive and nonreactive minerals and found a nonlinear impact of mineral composition on rock weathering rates. However different patterns of discontinuities in rocks were not examined. Here, we chose to simulate the weathering of mono-mineralic rocks with different types of discontinuities. The model simulates a far from chemical equilibrium system and does not consider any transport limitations.
These weathering conditions are representative of rocks and saprolite on or near the surface in non-arid environments.
The discontinuities in the model are intended to represent grain boundaries, joints, cracks, or stylolites that are partially filled with Simulation of fluid-rock interaction. The bulk rock components are marked in yellow, rock discontinuities in black, and fluid in white. Cross sections of the rock are shown at three stages of the simulation: (a) initial state, (b) step 220, and (c) step 221. Chemical weathering dissolves the rock minerals slower than the discontinuities between rock clusters. When a cluster is surrounded by fluid, it detaches from the surface and is removed from the simulation. Note that the black discontinuities dissolve more rapidly than the bulk rock [Colour figure can be viewed at wileyonlinelibrary.com] cement. Thus, the discontinuities have an intrinsic strength that binds the rock together, but they also dissolve more rapidly than the bulk rock, and this effect is included in the model.
The data from every simulation were saved as an object comprising all the information from the simulation, including the rock's initial properties and the dynamic properties of the rock in every step. These properties include an image of the rock in every step, a list of pixels that were dissolved, and the location and dimensions of detached fragments in every step. Using this object-oriented approach in MATLAB™, each simulation takes several minutes on a standard PC, and the data are uploaded into a MySQL database facilitating analysis of the datasets.

| Patterns of discontinuities
In our model, we used two kinds of discontinuity patterns: synthetic and natural ( Figure 2). Four different synthetic patterns were tested: (i) regular grid jointing; (ii) brick wall jointing; (iii) hexagonal jointing, simulating columnar patterns common in basalts; (iv; see also Supplementary Video S1) Voronoi tessellation, representing a coarsegrained crystalline rock. Weathering was also simulated for four natural rock patterns, obtained by binarization of outcrop images:

| Model calculations
In the initial state of our simulations, we define a grain or block as a region bounded by discontinuities. The width of discontinuities is constant and varies between one and three pixels for all simulated patterns. We also define the discontinuity density as the proportion of discontinuity pixels in the domain. For natural rock patterns in our simulations, this varies in the range of 2-30%, while for synthetic patterns this varies from 7% to 40%. For natural patterns, different values of discontinuity density were obtained by cropping the large natural rock images in different locations into the desired size (12 different initial images were used for perpendicular stylolites, 12 for parallel stylolites, 3 for diagonal cracks, and 10 for orthogonal cracks). For each cropped pattern image, we ran 6-10 realizations.
In simulations using synthetic patterns, discontinuity density is where ϵ is the discontinuity density, D is the intrinsic diffusivity of the discontinuity network, while D eff is the effective diffusivity through the bulk rock. We used the Tau Factor MATLAB™ application (Cooper et al., 2016) to calculate the tortuosity based on our 2D images.

| RESULTS AND DISCUSSION
3.1 | Impact of discontinuity density on rock weathering rates For all the rock patterns we tested, we found weathering rates to increase as the discontinuity density increased (Figure 3). This result is not surprising since the dissolution rate along the discontinuities is more rapid than the dissolution rate of the bulk rock. Moreover, this is consistent with field and experimental observations of weathering rates in fractured rocks (Eppes & Keanini, 2017;Røyne et al., 2008).
Our simulations also show that the type of discontinuity pattern has a significant impact on weathering rates (Figure 3), particularly at discontinuity densities <25%. For the synthetic patterns, at any given value of discontinuity density, the rates followed the order: grid > honeycomb > Voronoi > brick wall. In natural rock patterns, the order was less clear, although weathering in orthogonal cracks was faster than in diagonal cracks, and weathering in perpendicular stylolites was faster than in parallel stylolites. In addition, synthetic patterns generally weathered faster than natural rock patterns. This is probably due to the irregular nature of natural patterns and their inherently lower connectivity.
At discontinuity densities >25%, the weathering rates of all the patterns begin to converge. This may be because, at low discontinuity den-F I G U R E 3 Weathering rate as a function of discontinuity density in synthetic (red) and natural (blue) patterns. Each of the synthetic patterns shows a linear increase in weathering rate with the density of discontinuities. In the natural rock patterns, the trend is less clear. At lower discontinuity density, the weathering rate exhibits a strong dependence on the pattern [Colour figure can be viewed at wileyonlinelibrary.com] sities, the tortuosity of the pathways and their low connectivity acts as a limiting factor. As the discontinuity density increases, connectivity is expected to increase, facilitating the advance of the weathering front.
Although the discontinuity density is a critical parameter in determining weathering rates, our results suggest that additional parameters related to the geometry of the patterns are also likely to impact the way rocks weather. Specifically, for patterns in which the pathways are highly tortuous and poorly connected, rates are expected to be slower. F I G U R E 6 Snapshots of three simulations with different offsets and tortuosities. Each image shows the time step directly preceding a grain detachment event. The initial grain size is μ 0 = 2,160 pixels. Note the decrease in size of detached grains and in the difference in patterns of penetration of the reactive fluid into the rock. In addition, for the lowest level of tortuosity, the reaction front advances much more rapidly than at higher levels of tortuosity [Colour figure can be viewed at wileyonlinelibrary.com] F I G U R E 7 Weathering rate versus tortuosity in simulations with offset from grid (0% offset) to brick-wall (50% offset) and back to grid (100% offset), in three decreasing initial grain sizes μ 0 = 2,160 pixels (black), μ 0 = 234 pixels (blue), and μ 0 = 108 pixels (red) [Colour figure can be viewed at wileyonlinelibrary.com]

| Impact of tortuosity on weathering rates
In the simulations of synthetic rocks, each pattern type showed a decrease in weathering rate with increasing tortuosity (Figure 4a).
Moreover, the rates grouped into two distinct trends: (i) grid and brick wall, and (ii) honeycomb and Voronoi. This is probably due to the similarity in the geometry of the patterns within each trend. By contrast, for natural rock patterns, there is no clear dependence of weathering rate on tortuosity for individual pattern types (Figure 4b). This could be related to the irregularity and anisotropy of discontinuities in natural patterns, which can cause patterns with identical tortuosities to behave differently. In addition, the widely varying discontinuity densities in the natural patterns could also mask the apparent impact of tortuosity.
Each synthetic rock pattern showed a distinct trend with respect to mechanical weathering. Specifically, for each type, the proportion of mechanical weathering relative to the overall weathering increased with increasing tortuosity (Figure 5a). This relationship is likely to be the result of the dependence of tortuosity on grain size in our simulations; as grain size increases, larger grains are detached, increasing the proportion of mechanical weathering. However, in natural rock patterns (Figure 5b), no clear relationship between tortuosity and the contribution of mechanical weathering to the overall rate was found.
Moreover, the contribution of mechanical weathering to the overall rate was smaller in natural patterns than in synthetic patterns: in stylolite simulations the contribution of mechanical weathering was between 9% and 30%, while in synthetic patterns mechanical weathering contributed 30-47% of the overall weathering. This could be due to the complex patterns in rocks containing stylolites, in which blocks may be poorly connected to each other, resulting in a higher level of chemical weathering for each particle before it detaches, and an overall lower contribution of grain detachment.
To isolate the impact of tortuosity, we conducted a numerical experiment with simulations of synthetic patterns in which tortuosity changed systematically while maintaining the same level of discontinuity density (Figure 6). Starting with a regular grid, we introduced an offset in alternating layers to create brick wall patterns, which increased the tortuosity. In this method, the tortuosity varied from 1.95-2.75. In each offset, we ran six simulations with three different initial grain sizes: 2,160, 234, and 108 pixels.
Our results show a near-linear decrease in the weathering rate as the tortuosity increases from 1.95 to 2.75 (Figure 7) for all three grain sizes tested. Overall, the reduction in rate was 33%, 21%, and 27% for the 2,160, 234, and 108 grain size simulations, respectively. This significant effect means that, in addition to mineralogy and grain size, (b) F I G U R E 8 Impact of discontinuity density on mean detachment size in (a) synthetic and (b) natural patterns. The synthetic patterns show a nonlinear reduction in detachment grain size with increasing discontinuity density. In the natural patterns, the discontinuity abundance is less systematic, but the trend is similar [Colour figure can be viewed at wileyonlinelibrary.com] Mean detachment size as a function tortuosity in offset simulations for three initial grain sizes: (a) 2,160 pixels, (b) 234 pixels, and (c) 108 pixels [Colour figure can be viewed at wileyonlinelibrary.com] the tortuosity of the discontinuity pattern is likely to be a critical factor in determining the weathering rate in real rocks.

| Impact of discontinuity density and tortuosity on size of detached grains
We found that the mean detached grain size decreases nonlinearly with increasing discontinuity density for both synthetic and natural patterns (Figure 8). For the synthetic patterns, the detachment grain size drops by approximately 90% as the discontinuity density increases from 8% to 25% (Figure 8a). For natural patterns, there is a significant level of variability and the trend is far less clear (Figure 8b).
This is most likely a result of the differences between the initial conditions in the synthetic patterns and those in the natural patterns: in the synthetic patterns, the initial grain sizes are similar for any given discontinuity density, while in natural patterns, the initial grain size varies significantly.
The overall reduction of the mean detached fragment size with increasing discontinuity density is caused by two factors. The first is that increasing discontinuity density leads to a reduction in the initial grain size, which results in smaller detached grains. The second is that, as the discontinuity density increases, the chemical weathering rate also increases, causing the grains to undergo more dissolution prior to detachment.
To test whether tortuosity plays a role in the size of detached grains, we analyzed the results of the offset experiment described in Section 3.2 and found that the mean detachment size decreases with increasing tortuosity (Figure 9). This is because, in patterns with higher tortuosity, chemical dissolution has a longer time to act and reduce the size of the grains prior to detachment. This effect can be seen in the simulation snapshots in Figure 6: detaching grains in the grid simulation are larger than the detaching grains in the offset simulations.

| Impact of textural patterns on the size distribution of weathered grains
In all the rock patterns we tested, the grain size distribution of  (Figure 10d). For the natural rock patterns that we tested (stylolites and cracks), the initial block size distributions were multimodal. However, the size distribution of the detached fragments showed reduced modality ( Figure 11). Our results are consistent with the findings of Palomares et al. (1993), who showed that rocks with similar initial grain sizes fragment mainly along their uniformly distributed discontinuities, thus providing grains of uniform size, in contrast to rocks with anisotropic fabrics that do not disintegrate uniformly.
Grain size distributions of weathered grains strongly influence soil permeability and soil erosion (Cohen et al., 2015). Soils with a wide range of grain sizes are less permeable and erode less readily than soils with uniform grain size distributions (Cohen et al., 2015). Thus, we expect rocks with initial grid-like, or honeycomb discontinuity patterns to produce relatively uniform grain size distributions that form soils with higher permeabilities. By contrast, rocks with stylolites and cracks might be expected to produce soils that form impermeable layers.
Although there is significant variability, the grain size distribution of many soils and sediments often has a log-normal distribution (Gardner, 1956;Wagner & Ding, 1994). In our simulations, the only pattern that weathered into fragments with log-normal distributions is the Voronoi pattern. These patterns are common in the polycrystalline rocks that provide much of the weathered material to sediments, and it is likely that the log-normal distribution in sediments is influenced by the initial grain size distribution of the weathered rock. However, transport processes also strongly affect the size distributions of sediments, (Hunt & Sahimi, 2017), and we therefore expect the discontinuity patterns to have the strongest impact on the distribution of sediments that are relatively close to the source rock, such as in fluvial fans.

| CONCLUSIONS
We used a numerical model that incorporates both chemical and mechanical weathering to investigate the impact of rock texture on the weathering rate and the size of detached grains. Our results indicate that the weathering rate increases with increasing densities of discontinuities in the rock. We also found that increasing the tortuosity of the patterns leads to decreasing weathering rates.
Moreover, we found a strong impact of texture on the detached grain size distribution, and that higher discontinuity densities lead to smaller detached blocks. This has practical implications for risk assessment near cliffs or stone edifices: rocks containing stylolites with spacings of several centimeters could present less of a risk than rocks containing fractures with spacings of tens of centimeters.
The model we present here is a preliminary attempt to simulate the combined effects of chemical and mechanical weathering, and we can identify some limitations to our approach. Our simulations compare textures of different scales: the individual grains comprising a rock are often micrometer or millimeter in scale, while joints and stylolites are often present at the centimeter and meter scale.
Moreover, the time and spatial scales in the model are at present arbitrary, which severely limits its predictive power. In the current model, we simulated rocks with discontinuities that undergo faster weathering than the bulk rock. Although this assumption may be good for many scenarios, it will not be appropriate when the fractures are filled with an unreactive cement, such as quartz.
Calibrating the model, however, requires reliable field data, which are difficult to obtain because of the long time scales associated with weathering. Moreover, the difficulty of measuring the tortuosity of rock patterns in the field, as well as determining chemical and mechanical weathering rates, suggests that experimental work may provide a valuable alternative approach. Future work that focuses on improving the model by comparison with experimental and field-based measurements could provide solutions to some of these challenges.