A Statistical Approach To Explain the Solution Stoichiometry Effect on

Article pubs.acs.org/crystal

A Statistical Approach To Explain the Solution Stoichiometry Effect on Crystal Growth Rates Helge Hellevang,*,†,‡ Beyene G. Haile,† and Rohaldin Miri† †

Department of Geosciences, University of Oslo, P.O. Box 1047, Blindern, NO0316 Oslo, Norway The University Centre in Svalbard (UNIS), Pb. 156, 9171 Longyearbyen, Norway

‡

ABSTRACT: The solution stoichiometry strongly affects growth rates of ionic compounds, but it has been demonstrated that classically derived analytical growth rate models can only reproduce some of the observed variations. These models do not take into account surface diffusion, and a new statistical model was therefore constructed and solved to see if mechanisms, such as surface transport and ion complexation, can explain some of the discrepancies between the classical models and observed data. We showed that all experimental data can be reproduced, and the statistical model suggest that surface transport and ion complexation are required to obtain the narrow Gaussian/Lamé curves as observed, for example, AgCl. This was explained by the easier surface transport of the neutral complex, and growth units (charged ions) may then be supplied to the growth sites when the complexes are close to surface growth sites. Further numerical or laboratory tests are required to verify if this complexationsurface transport model can indeed be the explanation for some of the variation observed in the experimental data.

■

Z2/ρ and thereby the dehydration frequencies of the cations. This was true for the Me2+ series (Ba2+, Ca2+, Mg2+), but the AgCl salt behaved very differently with a much sharper Gaussian at the same Z2/ρ. The reason for this is unclear, but the strongly covalent nature may play a role. It was furthermore found that the Gaussian bells are generally symmetric, suggesting that the symmetric dehydration of both the surface sites and growth units (of opposite charge) plays an important role. There are several examples of mechanistic analytical models to explain the growth rate dependencies on the solution stoichiometry.13−15,18 These models works well to explain experimental data of some ionic compounds, but they are not general and cannot predict the range of stoichiometry dependencies observed for many compounds.11 In this paper we aim at improving the understanding of the mechanisms being responsible for the rates of crystal growth through a statistical approach. In contrast to the other mechanistic analytical models, we explore to what extent ion mobility along the mineral surface may affect mineral growth rates and if these

INTRODUCTION Understanding the factors that control crystal growth kinetics is of large importance for the production of a range of materials, from amorphous to nanosized crystal assemblages, to large inorganic pure to superpure crystals for industrial use.1−3 Understanding growth kinetics is especially important for the synthesis of proteins or other materials with functional groups, where crystal size, shape, and number and nature of defects are kinetically controlled.4−7 Crystal growth kinetics from aqueous solutions depend on factors such as solution oversaturation, temperature, and solution composition, and several mechanisms may operate at the same time.8−10 One feature that has received increased attention in recent years is the effect of solution stoichiometry on crystal growth rates. At fixed aqueous oversaturations, growth rates of ionic compounds (e.g., CaCO3, AgCl, BaSO4, etc.) have been shown to depend on the relative amount of the constituent cations and anions in solution.11−16 Also the advancement rates of individual crystal growth orientations are strongly affected by the solution stoichiometry.14,17 Rates tend to peak at or close to the stoichiometric solution, but the reduction from the peak rates for nonstoichiometric solutions varies from compound to compound.11 Hellevang et al.11 found that the variation in shapes (width) of the resulting Gaussian rate curves apparently depends on the charge density © XXXX American Chemical Society

Received: October 14, 2015 Revised: January 26, 2016

A

DOI: 10.1021/acs.cgd.5b01466 Cryst. Growth Des. XXXX, XXX, XXX−XXX

Crystal Growth & Design

Article

Figure 1. Sketch of the conceptual (statistical) model (a) complexation step (b) growth step.

vector sets rA and rB from the vector set ψ which correspond to growth unites A and B, respectively:

mechanisms can explain the Gaussian rate curves for some ionic compounds. Surface-transport through diffusion was considered as a rate controlling mechanism in the BCF growth theory and later related models,19,20 but diffusion were considered too slow in some later models, and direct transport to the growth sites from the bulk solution and subsequent attachment to the surface growth sites were instead assumed to be rate controlling.13,15,18,21 Recent molecular-scale simulations done on BaSO4 and NaCl support reactions over several transition states and that the rate-controlling step for these compounds do not involve surface diffusion.22,23 Nevertheless, some recent experimental observations and mathematical/ mechanistic models still point toward a surface diffusion control on crystal growth,21,24−27 and the question is now more on the possible distances covered by ions along the surface (nanometer scales) and whether the surface diffusion is rate limiting. To explore the effect of surface transport we defined a crystal surface with a finite number of growth sites at random position, and with growth units A and B also being placed at random positions. Only A and B can be incorporated into the surface, and formation of AB complexes therefore reduces the number of A and B that participate in growth, and growth rates are correspondingly reduced. Whether complexation occurs or a unit is taken up on a surface growth site is determined by the distance between the growth units and surface growth sites (a measure of the surface mobility), and the maximum distance for reactions to occur is varied. The model is run repeatedly until a statistical average of the growth rate is obtained. This statistical approach is shown to be capable of reproducing the dependencies of solution stoichiometry of all recorded ionic compounds (CaCO3, AgCl, MgC2O4xH2O, CaC2O4xH2O, CaSO4·x2H2O, BaSO4), and the individual contributions of aqueous complexation, growth site density, and surface mobility on the shape of the Gaussian bell is demonstrated.

rA = {ri|ri ∈ ψ , i ∈ Ii} rB = {rj|rj ∈ ψ , j ∈ Ij}

(1)

where Ii = {1, ..., NA} and Ij = {1, ..., NB} are index sets and rA∩rB = ϕ. The proposed statistical model is composed of two steps: (1) ion complexation and (2) surface uptake of growth units (mineral growth). These two steps are illustrated in the conceptual sketches in Figure 1. On the basis of the complexation step, only certain molecules in the random sets (i.e., rA and rB) would have the probability of the reaction based on the relative distance between unlike pairs; i.e., the distance between A and B should be less than a critical reaction radius, σreact. Considering this criterion, reactive molecules were sampled out as a new subset such that rreact = {ri ∈ rA , rj ∈ rB|drri j < σreact}

(2)

where drirj is the distance between unlike pairs in ψ: d(ri , rj) =

(xi − xj)2 + (yi − yj )2 where i ∈ Ii: j ∈ Ij (3)

Furthermore, among all members of rreacthereafter referred to as a populationa subset of individuals growth units are chosen based on a simple random sampling (SRC) technique; i.e., each growth unit is chosen randomly and entirely by chance via a probability (complexation probability PC) which remains constant during the sampling process. This has been done to take the formed complexes out of population and make a new population rp with individuals capable of incorporating into the growth sites. rp = {rn|P(rn ∈ rreact) = PC , n ∈ In}

■

(4)

where In = {1, ..., Nreact} and Nreact is the number of growth units, A and B, that can potentially be taken up by the surface, given that they are sufficient close to a surface reaction site. Once a set of molecules with potentials of growth is extracted from the complexation step, the algorithm checks if growth units A or B can be incorporated into the surface. Growth surface is composed of two different growth sites termed GA

GROWTH RATE MODEL We consider a two-dimensional system composed of N = NA + NB growth units of type A and B interacting with a crystal surface of size Ns. We describe all possible coordinates for growth units with a random vector set ψ⊂R2. Two random B



Article

Figure 2. (a−c) Variation of the mineral growth rate (normalized rate f ′) as a function of the solution stoichiometry (ξ) independently varying the reaction radius (Θ) and the surface density of growth sites (Λ), at a fixed probability of ion complexation of zero. All parameters are normalized and in dimensionless form.

and GB. Let us choose two random vector sets rGA and rGB from the sample space ψ to define the position of GA and GB, respectively:

to size of the system (see nomenclature). Essential steps of the model are summarized in the following algorithm:

rGA = {rk|rk ∈ ψ , k ∈ Ik} rGB = {rm|rm ∈ ψ , m ∈ Im}

(5)

where Ik = {1, ..., NGA} and Im = {1, ..., NGB} are index sets and rGA∩rGB = ϕ. This reaction step implies that only certain molecules in the random sets (i.e., rp) would have probability of incorporation into the surface based on the relative distance between ionic growth units and growth site of the opposite kind; i.e., the distance between A and GB should be less than a critical growth radius, σG and similar for B and GA. These rules summarized as 1. Growth units and growth site must be of the same type 2. Growth units must be located in the growth radius of growth site Considering these two rules, those molecules that incorporate into site GA and GB were sampled out from rp as a new subset rG: ⎧ < σG ⇒ rl ∈ rG , l ∈ Ip ⎪ rl ∈ rp ∩ rA , drr l GA ⎨ ⎪ ⎩ rl ∈ rp ∩ rA , drrl GB < σG ⇒ rl ∈ rG , l ∈ Ip

(6)

Because the shape of some of the experimental data could not be well fitted by Gaussian curves, we instead used the Lamé (superellipse) function to compare modeled and experimental data:28

where Ip = {1, ..., Np} and Np is size of the vector set rp. We call rG as the growth set, and the size of this set gives the number of incorporated molecules g. This process should be repeated for sufficiently large n steps, while all the model parameters remain unchanged, with the exception of the locations of growth units which are randomly updated. Growth rate is defined as an ensemble average such as G = ⟨g ⟩ =

1 n

log ξ a

i=1

+

f′ b

n

=1

(8)

where ξ is the concentration ratio of cation and anion in solution, and f ′ is the normalized growth rate or growth frequency. Doing this, we can easy quantify the shape of the model and experimental curves using the value of the Lamé exponent n. To fit Lamé curves to experimental data, we fixed a = 1 and b = 1 and obtained n for the different data sets using a Trust-Region quasi-Newton algorithm. An uncertainty of the

n

∑ gi

n

(7)

In order to keep the statistical model independent of the size of system, the model parameters were normalized with respect C



Article

Figure 3. (a−c) Variation of the mineral growth rate (normalized rate f ′) as a function of the solution stoichiometry (ξ) independently varying the probability of ion complexation (Pc) and the surface density of growth sites (Λ) at a fixed reaction radius (Θ) of 0.3. All parameters are normalized and in dimensionless form.

obtained n-values of ±10% was used for all minerals, except for calcite where several experimental data series are available and the ranges in n-values from these were used instead.

growth sites (Figure 3). The sensitivity of the peaks to the fraction of surface growth sites increases at Pc values lower than 0.5 and reaches a maximum at Pc = 0.0. The general pattern of a narrowing of the peaks at higher amounts of ion complexation comes from an increasing difference between the concentrations of the surplus and deficient ions in solutions; i.e., the concentration ratio ξ of the charged ions close to the mineral surface will deviate more from unity than the corresponding total concentration of these. To illustrate this we can use a simple example where we have 100 A and 10 B ions, giving a A/B ratio of 10. If five A ions complexes with B, the true ratio of the A+ and B− ions that can potentially be taken up at the surface has changed to 95/5 = 19, further away from the stoichiometric 1:1 ratio. The simulations up to now generated symmetric peaks as the same reaction radiuses were defined for both growth units, whereas experimental data for some ionic compounds and at some conditions suggest asymmetric curves and peak rates at nonstoichiometric solutions (Hellevang et al., 2014). In the model, only the reaction radiuses for the growth units can create asymmetries, and we therefore varied Θ individually for A and B (Figure 4). In this modeling we kept Λ constant at 10−4, ΘA constant at 0.2, varied ΘB from 0.005 to 0.2, and two sets were run at Pc of 0.333 (Figure 4a) and 0.0 (Figure 4b). The simulations show that ion complexation reduces the degree of asymmetry. Without ion complexation (Pc = 0.0), rates peak at an ionic ratio A/B of less than 0.01 at a ΘB/ΘA ratio of 0.025 (Figure 4a). At Pc = 0.333, rates peaks at an ionic ratio of 0.1 at the same ΘB/ΘA value (Figure 4b), and finally at a Pc of 1.0 rates peak at a A/B ratio of 0.4 (not shown). The reason for the reduced asymmetries with increasing Pc is that the number of A and B being removed from the growth assemblage by ion complexation is determined by the largest of the Θi, and a large part of the ions with low Θ will therefore be removed prior to the surface growth step. Limiting the Parameter Space of Ionic Compounds. Experimental data covering a range of ionic compounds were reviewed in Hellevang et al.,11 and the shape of the Gaussian curves were analyzed using an empirical model. It was shown that AgCl16 and BaSO429 showed the sharpest and widest Gaussian bells respectively, and all the other reviewed ionic compounds (CaCO3,12,17,30 MgC2O4xH2O,31 CaC2O4xH2O31 CaSO4·x2H2O32) had curves that were intermediate between these two. To compare the shapes of the experimental data with

■

RESULTS AND DISCUSSION Modeled Growth Curve Variations. The variation in growth rates with ionic activity ratio r = A:B was simulated for a system with ionic reaction product Ω = 9 × 10−4 and varying Θ (reaction radius), Λ (surface fraction of growth sites) and Pc (probability of ion complexation). Rates at different parameter sets varied, and all rate data were therefore normalized by dividing by the maximum growth value for each set at the stoichiometric aqueous solution. The normalized rates (growth frequencies) are reported as f ′. The first simulations were done without taking into account ion complexation (Pc = 0) and varying Λ and Θ (Figure 2). The same values of Θ were used for both types of growth units. Using a low fraction of surface growth sites (Λ = 10−6) leads to curves that go from a classical Gaussian at the low Θ value (0.1) to a broader bell-shaped curve at larger values of Θ (Figure 2a). The broad bell-shaped curves come from the relatively large number of growth units compared to the growth sites. With large reaction radiuses, more molecules will be available per growth site, and the stoichiometry dependency is reduced. The stronger stoichiometry dependency at smaller reaction radiuses comes from an increasing local environment dependency on the growth (the deficient growth component can only attach to the local growth sites and not choose from the global population of growth sites). Increasing the fraction of surface sites to 10−4 results in a general narrower and sharper peak, and reduces the variation of the curves with Θ (Figure 2b). At a further three times increase in the fraction of growth sites, rates become independent of Θ (Figure 2c). The curves are strictly speaking not Gaussian, perhaps only for the combination of few surface sites and small reaction radius (Figure 2a), and curves can be better represented by Lamé curves. The second set of simulations was done to see the effect of ion complexation. In these simulations the reaction radius was kept constant at 0.3 for both growth units, and the fraction of growth surface sites was varied from 10−6 to 10−4 (Figure 3). Simulations with Pc = 0.0 correspond to the ones shown in Figure 2. There is a general shift in the width of the peaks with increasing ion complexation. At Pc of 0.5 or higher, peaks are very narrow and do not vary much with surface density of D



Article

Figure 4. (a, b) Variation of the mineral growth rate (normalized rate f ′) as a function of the solution stoichiometry (ξ) independently varying the probability of ion complexation (Pc) and the ratio of the reaction radiuses of the cation (ΘA) and anion (ΘB). All parameters are normalized and in dimensionless form.

the numerical curves found in this study (Figures 2−4), we used the Lamé function (eq 8) and the Lamé-coefficient n obtained by regression for experimental growth rate curves for solution stoichiometries (log ξ) in the range −1.0 to 1.0. For AgCl and BaSO4 we found n to be 0.445 (±0.014) and 1.73 (±0.058) respectively (Figure 5). These values, representing the whole range in experimental shapes, were used as lower and upper boundaries for the comparison to the numerical model. We could then immediately filter out those simulations (model parameter sets) giving values outside this range, and this is shown in Figure 6. In this figure, the lame parameter n was calculated for discrete pairs of model parameters reaction radius (Θ), probability of ion complexation (Pc), and fraction of surface growth sites (Λ) and values of n in the interval 0.44 to 1.74 are illustrated with colors, whereas parameter pairs providing n-values outside the range is shown without color (black squares). Two general features can be observed from these plots: First, lower values of n is only found at larger reaction radiuses (Θ), and second, the variation in n-values are getting progressively larger with Θ (Figure 6). At the lowest reaction radius (Θ = 0.01), the Lamé n values are limited to the range from 1.13 to 1.46, and there are many possible combinations of Pc or Λ for a minerals such as calcite and CaC2O4xH2O (shown below). Parameter sets (Pc, Λ, Θ) for the minerals with the largest or smallest n-values (such as AgCl and BaSO4) may be better constrained than those with n values in the medium range. It is also easy to see that high Θ values (0.3 and 0.5), representing large surface mobility, have n-values for some of the parameter sets that are outside that found for the entire range in experimental data. For example, at Θ = 0.5, if the probability for ion complexation Pc is 1.0, n-values between 0.44 and 1.73 can only be found at the lowest number of surface sites (Λ = 10−6). Figure 6 indicates that different minerals have specific ranges in their possible parameter sets. This can be further visualized by printing the possible parameter sets for each individual mineral. This was done by first estimating n-values and their statistical uncertainty range for all experimental growth rate data, and then filling in colors for those parameter sets that are within the n-values found for the mineral (Figure 7). Blank (no color) indicates parameter sets outside the range in n found for the mineral. If the n-values of the minerals do not correspond to any parameter set, but instead fall between two

Figure 5. Lamé parameter n obtained by regression from experimental data of BaSO452 and AgCl,16 providing the total range of n-values for all experimental data (these two minerals represented the end-member shapes).

of the parameter values, colored lines are indicated for that border (this is mostly the case for the sharp changes in n at the higher Θ, e.g., for BaSO4 (Figure 7a) and MgC2O4xH2O (Figure 7c) at Θ = 0.5). It is clear from this plot that some minerals are only found within a quite narrow parameter space (e.g., BaSO4, AgCl, MgC2O4xH2O), whereas other minerals (such as all the Ca-minerals) have growth rates that can be explained by a very large range in parameter values (Figure 7d−f). AgCl has been shown to have a very sharp drop in rates as the solution changes from the stoichiometric 1:1 Ag+:Cl− ratio.11,16 This translates into a very low n-value (0.34 ± 0.034), only being found at the highest Θ value of 0.5 (Figure 7b). Furthermore, these low values can only be found at this Θ value if the probability of ion complexation is 0.33 or higher. The relatively low n of MgC2O4xH2O (0.691 ± 0.069) similarly constraints the parameter space to the larger values of Θ and also here the ion complexation is indicated, but the total number of possible parameter sets is larger (Figure 7c). The Ca minerals generally cover much wider parameter spaces than suggested for the Ba, Ag, and Mg compounds (Figure 7d−f). This is especially true for calcite. If we look at the total number of variation in n (0.776−1.505), any combination of Pc and Λ at E



Article

Figure 6. A map of the variation of estimated Lamé n-values for different discrete pairs of model parameters (reaction radius (Θ), probability of ion complexation (Pc), and fraction of surface growth sites (Λ)). Black (no colored) squares represent modeled values outside the range found for the experimental data.

the lower values of Θ (0.01 and 0.02) will give n-values in this range (Figure 7f). There are also a good number of possibilities for Θ up to 0.1. There is however some variations in the type of data used. Data from bulk-rate experiments cluster in the middle, but are also present at the higher Θ values, whereas the AFM single-face experiments are more spread but mostly found in the lower Θ values (Figure 7f). Capabilities of the Kinetic Ionic Ratio (KIR) Models To Predict the Range in Observed Growth Rate Data. Peak growth rates at close to stoichiometric solutions have earlier been explained in several analytical mechanistic models.13−15,18 The main framework of these models was developed in the “kinetic ionic ratio” (KIR) model by Zhang and Nancollas,15 and based on the assumption that the rate-controlling steps were surface attachment/detachment of ions, whereas solute transport to the surface sites was treated as fast relative to the surface processes and thereby not rate limiting. As a basis for the model the following assumptions were made: (1) Attachment is proportional to the ion activities in solution, i.e., αA = kA[A] and αB = kB[B], where α denotes attachment frequency, k is a first-order rate constant, and [] denotes activity. Attachment rates are identical for all kink site. (2) Detachment frequencies are related to the binding energy of surface ions, and they do not depend on the solution ion activitites. (3) Detachment frequencies are the same for both ions, since the same bonds are broken (vA = vB = v). (4) Detachment from non-kink sites is negligible for supersaturated aqueous solutions. Furthermore, the solubility product KS should equal the ratio between detachment and attachment frequencies (KS = vAvB/kAkB), and growth rates can be predicted analytically only

knowing the kinetic ionic ratio (ξ′ = kA[A]/kb[B])15 and the saturation index. The background for this is the assumption that the rate is limited by one single and reversible elementary reaction step. This assumption has however recently been disproved for minerals such as BaSO422 and NaCl,23 and cannot be regarded as generally true for all minerals. There have been some discussions on the theoretical background for the attachment frequency of ions,33 but in the mechanistic and analytical models discussed here, dehydration of cations are assumed to control the growth rates. Dehydration frequencies of anions are much larger, implying that, to get the symmetric Gaussians as often observed in the experimental data, dehydration of an ionic growth unit A (cation) must be followed by a corresponding dehydration of a surface growth site B (anion), and vice versa.13,14,34 This may however be further complicated by indications that the order of the breaking of bonds in the detachment steps is important for the rate.22 In recent years, access to new experimental data has resulted in new insight in the KIR model, and also in modification. Stack and Grantham18 modified the equation for nucleation rates and used a Newton−Raphson minimization technique to fit the model parameters (attachment and detachment frequencies) to experimental AFM data on calcite growth. Their work suggested that attachment and detachment frequencies are not generally the same for Ca2+ and CO32−. In support of their model they found that their estimated ratio of averaged detachment and attachment frequencies was in good agreement with the solubility product of calcite. Larsen et al.14 compared modeling results using KIR to experimental data from AFM and found that only the obtuse step velocities could be found using realistic and fitted model parameters. No fit could be found for the acute step velocities, giving a narrow Gaussian and with peak rates deviation from the stoichiometric 1:1 F



Article

Figure 7. Possible parameter spaces (Pc, Λ, Θ) for the different minerals. Colored squares show discrete parameter sets with estimated Lamé n-values within the range estimated for the mineral.

zCa2+/CO32− ratio. The skewed Gaussian points to the difficulty of incorporating the large CO32− ion compared to

Ca2+ in the structurally closed kink sites at the acute steps. This observation suggests that the assumption in the KIR model of G



Article

Figure 8. Gaussian curves estimated with Zhang and Nancollas (KIR) model15 for (a) step movement/spiral growth; and (b) diffusion-controlled growth, and varying the saturation ratio S from 1 (equilibrium) to 10 (100 times supersaturation). Translated into the Lamé n-values, these shapes correspond to n = 1.904 ± 0.044 (surface reaction, S = 2) and n = 1.097 ± 0.053 (transport controlled, S = 2).

Figure 8. The first we notice is that the Gaussian (correction) for transport-controlled growth (Figure 8b) is much narrower than for surface controlled growth (Figure 8a). We also notice that the variation is very limited within each of these growth domains. Note that diffusion-controlled growth is likely to only occur at low aqueous supersaturations,8,19 and the highest value of the saturation ratio is therefore not very relevant for this process. Translated into the Lamé n-values, these shapes correspond to n = 1.904 ± 0.044 (surface reaction, S = 2) and n = 1.097 ± 0.053 (transport controlled, S = 2). We see that the value for surface reaction is larger than any of the estimated values from experimental data (Figure 6, 7), whereas for the diffusional transport is in the range for all the Ca-minerals (Figure 7d−f). The values are higher than estimated for AgCl (Figure 7b) and MgC2O4xH2O (Figure 7c), implying that the drop in rates at nonstoichiometric solutions are much sharper for these two minerals than possible for the original KIR model. Some of the basic assumptions of the KIR model were changed in the Stack and Grantham model,18 and optimal values of attachment and detachment rates were found by regression, and detachment/attachment rates can therefore be varied independently. The range of the KIR model results can therefore be expanded further (Figure 9). By doing so, the model can predict any (narrow) shape, also the Gaussian observed for AgCl. It would be of interest to perform similar molecular scale simulations as done in Stack et al.22 to check if the regressed values of the Stack and Grantham model compare well also with detachment/attachment frequencies of ionic compounds such as AgCl. Also the Wolthers et al. model13 produces broad Gaussians, similar to the two former, and the freedom of varying the shape is constrained by the solubility product. The only difference is that their complexation model provides a pH-dependency of the positions of peak rates. The Wolthers et al. model13 can also therefore not reproduce the narrow Gaussians observed for minerals such as AgCl and MgC2O4xH2O. One additional challenge with the model is that one of the detachment frequencies are used as a fitting parameter, implying that the attachment and detachment frequencies lose their required dependencies. Adsorption or Surface Diffusion Control? Numerous growth rate studies have observed parabolic rate laws of the form v = k(S − 1)2 for mineral growth done at low supersaturations.8,36−38 The background for parabolic rate

identical ion attachment rates at any kink site does not hold, and a more complex model taking into account variability of kink geometries must be developed.14 Wolthers et al. improved the KIR model further two years later and added a surface complexation model for calcite, taking into account surface speciation and the pH dependency of attachment and detachment frequencies of growth units. The surface complexation model is an important improvement to the KIR model as the changes in surface charge have been proposed to affect surface diffusion and incorporation rates of the cations and anions.35 In the model, growth rates of CaCO3 peak close to the stoichiometric ratio at pH 10.15, whereas peak rates are shifted far into CO32− deficient solutions at lower pHs. This was explained by an increasing reliance of growth on bicarbonate rather than carbonate. A comparison between the model results and various experimental growth rate data at different pHs has later left some doubt on the validity of the model, but it was not possible to draw any conclusions as no systematic data set exists for a proper model verification.11 Most verifications of the KIR model have been done using growth rate data from calcite. Larsen et al. showed that the model (the original one from Zhang and Nancollas) can reproduce some but not all of their AFM growth rate data.14 Also the process-based model of Wolthers et al.13 could reproduce some of the data, but there is a very large scatter when compared to the “global” calcite data set (see their Figure 7a). It is interesting to drive the comparison further and see if the KIR models can reproduce the range in Gaussian shapes for all phases reviewed in Hellevang et al.11 (CaCO3, AgCl, MgC2O4xH2O, CaC2O4xH2O, CaSO4·x2H2O, BaSO4) and also illustrated in Figure 7. We start out checking the degree of possible variability in the two simplest models, the original from Zhang and Nancollas,15 and the one from Stack and Grantham.18 In the original model, attachment and detachment frequencies are dependent through the ion solubility product, and only the kinetic ion ratio and saturation index can be adjusted.15 The correction functions for two phenomena were shown by Zhang and Nancollas, normal step movement/spiral growth being surface controlled and transport controlled growth (bulk solution diffusion → adsorption → surface diffusion). The correction functions at saturation ratios (here defined as S = ([A]α[B]β/Ks)1/2) from 1.0 (equilibrium) to 10 (100 times supersaturated) is shown in H



Article

expected that the charge, water density, and water selfdiffusivities will show a wave-like distribution with distance from charged mineral surfaces.49 This will imply that cations and anions should also be similarly distributed in zones (according to charge) out from the mineral surface, and it would therefore be surprising if the cations and anions have very similar lateral surface mobility over charged surfaces (this is required if surface diffusion is rate controlling and the effect of solution stoichiometry on rates, the Gaussian, is symmetric). Nevertheless, it is clear that both theoretical models and experimental observations have been used to support surface diffusion as an important control on mineral growth rate at low supersaturations. Suggested “New” Mechanisms To Explain the Effect of Solution Stoichiometry on Mineral Growth Rates. The background for making a statistical model was to see if transport and ion-complexation may influence the effect of solution stoichiometry on growth rates of ionic compounds. From the modeling work we found strong links between the falloff of growth rates from their peak rates, and parameters such as probability of ion complexation (Pc), the reaction radius (Θ), and the number density of growth sites (Λ). The narrower curves (small Lamé n-values) of minerals such as AgCl and MgC2O4xH2O could only be found at large reaction radiuses in combination with ion complexation (Figure 7b,c). These mechanisms are ignored in the kinetic ionic ratio (KIR) models as they assume rapid transport to the growth sites (no transport limitation) and that the ions are present at kink sites at concentrations governed by the bulk-solution solubility product. Lamé or Gaussian curves that are symmetric around the peak rates imply that the mechanism responsible for growth is also symmetric; i.e., there is no difference in the incorporation rates of the cation at anion-deficient solutions, and vice versa.11 In our model, growth units can be transported by surface diffusion and being taken up by a growth sites (kinks) if they are within a certain distance Θ from the reaction site. Using different reaction radiuses for the cation and anion resulted in shifts of the peak rates to either cation or aniondeficient solutions and asymmetric curves, and the magnitude of the shift was also connected to the amount of ion complexation (Pc) (Figure 4). Although some of the mineral experiments may suggest such asymmetric and shifted curves, most, and especially the bulk-rate experiments, show highly symmetric curves close to the stoichiometric ratio of the mineral.11 It is therefore questionable if surface diffusion alone can explain the very narrow and highly symmetric Lamé/ Gaussian curves of minerals such as AgCl.11 Instead, as supported in the present modeling, symmetric narrow Lamé/ Gaussian curves could hypothetically result from ion complexation, affecting the transport properties of the ions along step terraces. Ion complexation results in the formation of neutral species Aα+ + α/βBβ‑ = ABα/β, and both cations and anions can be transported in the complexes with the same surface diffusivity. These complexes are also likely transported more easily (and longer) along the surface as they are neutral and will be less affected by the mineral surface charges/forces. As complexes are transported they will decompose and reform according to a thermodynamic equilibrium, dictated by the forward/backward kinetics of the elementary reaction. This will ensure that the individual ions, A and B, will be available for uptake when the complex is decomposed sufficiently close to the growth sites. This hypothesis is further supported by large electronegativity of Ag+, being much larger than, e.g., Ca2+, giving the highly covalent nature of the ion.50 This is reflected

Figure 9. Gaussian curves estimated with the kink-site nucleation and propagation model of Stack and Grantham,18 using their estimated parameters, but using a range in detachment/attachment frequency ratios.

laws is that the rate-determining mechanism is spiral growth.8,19,33 Whether rates are controlled by surface diffusion or by direct adsorption is determined by some fundamental properties of the mineral surfaces and surrounding solution. Because of the relatively large total surface area compared to the area governed by kink sites, the majority of ions will be adsorbed at some distance from the kink sites.19,33 Whether the growth is diffusion controlled or adsorption controlled is determined by the average distance between growth steps and growth units. If the distance between growth steps y0 are on average (much) larger than the distance of adsorbed molecules from kink (growth) sites xad, then the growth will be diffusion controlled. If, on the other hand, y0 ≪ xad, the growth will be adsorption controlled.33 Because surfaces on molecular scales are not uniform and dynamic, there will also be cases where a mixture of the two end-member mechanisms is controlling the growth. On the basis of a review of a number of electrolytes growing according to parabolic rate laws, Nielsen concluded that the rate-determining step of these were governed by two simultaneous phenomena: dehydration of cations followed by diffusion of the dehydration cations into the growth position.33 There have also been several experimental studies pointing toward surface diffusion control on growth rates at low supersaturations.8,24,25,39 To further understand the mechanisms that control crystal growth, Kwon et al.21 developed a numerical model based on the Finite-Element Method and discussed their numerical observations in relation to the classical mechanistic and analytical models from Chernov,40 Gilmer, Ghez, and Cabrera,41,42 and van der Erden,43 and more recent approaches to crystal growth.44−48 From a comparison to Chernov40 (direct transport of growth units to growth steps, no surface diffusion), they concluded that it is impossible to separate the effects of pure attachment effects versus pure transport (bulk diffusion) on crystal growth rates. Their numerical model also agreed with the earlier conclusions that the relative importance of surface diffusion and adsorption (bulk diffusion) on growth rates strongly depends on the distance between growth steps (L), with small values of L leading to largely overlapping surface diffusion fields and reduced uptake by this mechanism. One argument against a surface diffusion control is that the mobility may be very low for adsorbed ions. This may be especially true for ions adsorbed in the inner high-water-density region close to the surface. It is I



Article

model (Figures 6 and 7) and that the parameter space was quite well constrained for some of the minerals. This was especially true for AgCl and MgC2O4xH2O, having the narrowest Gaussian/Lamé curves. The modeling suggested that these narrow curves, especially for AgCl, depend on a combination of large mobility of growth units and ion complexation. We hypothesized that individual charged ions may be prevented from being transported laterally because of charged surfaces, and that cations and anions are expected to show different surface diffusivities according to their affinities for the surface. Therefore, symmetric Gaussian/Lamé curves can only be explained by a mechanism where ions are complexed and transported, and complexes being in equilibrium with their ions leading to a supply of growth units when complexes are sufficient close to the growth sites. There are at present no models that can verify this hypothesis, but it is clear that laboratory experiments using, e.g., ligands that affect the ion complexation, should impact the Gaussian/Lamé curves.

in the significantly larger fraction of neutral complexes of AlCl in solution than e.g., CaCO30, CaSO4, and BaSO40 (Table 1). There is however only a modest correlation between the lamé parameter n and the complexation of ions into neutral species, and complexation values in Table 1 is shown for the bulk Table 1. Comparison of Lamé Parameter n and Ion Complexation (Neutral Complexes) for Compounds Where Thermodynamic Data Were Available (phreeqc v3 databases)51a aqueous complex 0

AgCl CaCO30 BaSO40 CaSO40

pH 7

pH 10

lamé n

0

AB /A

AB0/A

0.45 0.78 1.58 1.06

0.94 0.00076 0.24 0.11

0.94 0.38 0.24 0.11

■

a

Values of the concentration ratio of the complex AB0 over the cation A is shown for pH 7 and 10.

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected].

solution rather than close to the mineral surface in the diffuse boundary layer, and more work is therefore needed before we can arrive at conclusions about the relation between ion complexation and the transport properties. Growth rates of MgC2O4xH2O is also constrained to Θ values ≥0.1 and could also to some extent be dependent upon ion complexation and transport. Our simulations cannot constrain growth parameters for calcite well, as there are many parameter combinations that can explain the Lamé/ Gaussian curves (Figure 7f). The model also suggest other parameter relations, such as for BaSO4 where the Lamé curves for Θ values ≥ 0.1 can only be explained if the number of growth sites (surface fraction of growth sites Λ) is low. Such links are however very uncertain, and no conclusions will be attempted here. Another weakness of the statistical model may be that there is no direction dependency on the growth; that is, the observed differences in rates and mechanisms between the acute and obtuse steps14 cannot be resolved. Any directioncontrol on surface transport (surface diffusion) is however not expected, at least not for the neutral ionic complexes. The present statistical model result cannot be verified since there are no numerical or analytical models that can, at the same time, handle ion complexation, variation in surface diffusivity of complexed and noncomplexed ions (at different distances from the surface), and the classical step-velocity mechanisms. More work is therefore needed to understand if, and if so, to what extent, ion complexation and surface diffusion are mechanisms that actually impact the solution stoichiometry dependency of growth rates.

Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS We greatly appreciate the constructive reviews given by the editor and two anonymous reviewers. This work has been partially funded by the FME SUCCESS centre for CO2 storage under Grant 193825/S60 from Research Council of Norway (RCN). FME SUCCESS is a consortium with partners from industry and science, hosted by Christian Michelsen Research AS.

■

NOMENCLATURE ξ bulk solution activity ratio of cation/anion rA vector set holds coordinates of growth unit A rB vector set holds coordinates of growth unit B NA number of growth unit A NB number of growth unit B Nreact total number of reactive growth units (A + B) rp vector set holds coordinates of growth units after complexation step Np total number of reactive growth units after complexation step rGA vector set holds coordinates of growth surface sites of type A rGB vector set holds coordinates of growth surface sites of type B NGA number of growth surface sites of type A NGB number of growth surface sites of type B rG vector set holds coordinates of growth units incorporated in the surface g number of incorporated growth units = size of rG Θ reaction radius (dimensionless) Pc probability of complexation 0 ≤ Pc ≤ 1, providing an equilibrium A + B ↔ AB0 concentration G ensemble average of incorporated growth units (A + B) after n steps Ns total number of lattice units r coordinates (in lattice units) where rdi = Ni : rD

■

CONCLUSIONS Experimental data have earlier been reviewed and showed large variations in the width of the Gaussian (Lamé) curves for different minerals.11 It was also suggested that the classical analytical models cannot be used to explain the observed variations.11 These classical models assume that rate limiting step for growth is ion and surface dehydration, and they do not take into account surface diffusion. A new statistical model was therefore constructed and solved to see if surface transport, ion complexation, and surface density of growth sites affect the dependency of solution stoichiometry on mineral growth rates. We showed that all curves can be reproduced by this statistical

s

Θ= J

σ Ns

rD = {rdi|rdi ∈ ψ , i ∈ Ii} normalized reaction or growth radius DOI: 10.1021/acs.cgd.5b01466 Cryst. Growth Des. XXXX, XXX, XXX−XXX

Crystal Growth & Design Ω= Λ= Ξ=

■

NANB Ns

reaction product of growth units

NG Ns Ω Λ

surface fraction of growth sites

Article

(21) Kwon, Y.-I.; Dai, B.; Derby, J. J. Assessing the dynamics of liquid-phase solution growth via step growth models: From BCF to FEM. Prog. Cryst. Growth Charact. Mater. 2007, 53, 167−206. (22) Stack, A. G.; Raiteri, P.; Gale, J. D. Accurate Rates of the Complex Mechanisms for Growth and Dissolution of Minerals Using a Combination of Rare-Event Theories. J. Am. Chem. Soc. 2012, 134, 11−14. (23) Liu, L.-M.; Laio, A.; Michaelides, A. Initial stages of salt crystal dissolution determined with ab initio molecular dynamics. Phys. Chem. Chem. Phys. 2011, 13, 13162−13166. (24) Sleutel, M.; Vanhee, C.; Van de Weerdt, C.; Decanniere, K.; Maes, D.; Wyns, L.; Willaert, R. The Role of Surface Diffusion in the Growth Mechanism of Triosephosphate Isomerase Crystals. Cryst. Growth Des. 2008, 8, 1173−1180. (25) Sazaki, G.; Okada, M.; Matsui, T.; Watanabe, T.; Higuchi, H.; Tsukamoto, K.; Nakajima, K. Single-Molecule Visualization of Diffusion at the Solution−Crystal Interface. Cryst. Growth Des. 2008, 8, 2024−2031. (26) Kwon, Y.-I.; Derby, J. J. Modeling the coupled effects of interfacial and bulk phenomena during solution crystal growth. J. Cryst. Growth 2001, 230, 328−335. (27) Koutsopoulos, S. Kinetic Study on the Crystal Growth of Hydroxyapatite. Langmuir 2001, 17, 8092−8097. (28) Gielis, J. A generic geometric transformation that unifies a wide range of natural and abstract shapes. Am. J. Bot. 2003, 90, 333−338. (29) Kowacz, M.; Putnis, A. The effect of specific background electrolytes on water structure and solute hydration: Consequences for crystal dissolution and growth. Geochim. Cosmochim. Acta 2008, 72, 4476−4487. (30) Nehrke, G.; Reichart, G. J.; Van Cappellen, P.; Meile, C.; Bijma, J. Dependence of calcite growth rate and Sr partitioning on solution stoichiometry: Non-Kossel crystal growth. Geochim. Cosmochim. Acta 2007, 71, 2240−2249. (31) Chernov, A. A.; Petrova, E. V.; Rashkovich, L. N. Dependence of the CaOx and MgOx growth rate on solution stoichiometry. NonKossel crystal growth. J. Cryst. Growth 2006, 289, 245−254. (32) Zhang, J.; Nancollas, G. H. Influence of calcium/sulfate molar ratio on the growth rate of calcium sulfate dihydrate at constant supersaturation. J. Cryst. Growth 1992, 118, 287−294. (33) Nielsen, A. E. Electrolyte crystal growth mechanisms. J. Cryst. Growth 1984, 67, 289−310. (34) Christoffersen, J.; Christoffersen, M. R. Kinetics of spiral growth of calcite crystals and determination of the absolute rate constant. J. Cryst. Growth 1990, 100, 203−211. (35) Hong, M.; Teng, H. H. Implications of solution chemistry effects: Direction-specific restraints on the step kinetics of calcite growth. Geochim. Cosmochim. Acta 2014, 141, 228−239. (36) Larsen, K.; Bechgaard, K.; Stipp, S. L. S. Modelling spiral growth at dislocations and determination of critical step lengths from pyramid geometries on calcite {} surfaces. Geochim. Cosmochim. Acta 2010, 74, 558−567. (37) Boistelle, R.; Doussoulin, A. Spiral growth mechanisms of the (110) faces of octacosane crystals in solution. J. Cryst. Growth 1976, 33, 335−352. (38) Lundager Madsen, H. E.; Boistelle, R. Growth kinetics of the (001) faces of hexatriacontane (C36H74) in solution. J. Cryst. Growth 1979, 46, 681−690. (39) Land, T. A.; De Yoreo, J. J.; Lee, J. D. An in-situ AFM investigation of canavalin crystallization kinetics. Surf. Sci. 1997, 384, 136−155. (40) Chernov, A. A. The spiral growth of crystals. Soviet Physics Uspekhi 1961, 4 (1), 116. (41) Gilmer, G. H.; Ghez, R.; Cabrera, N. Analysis of combined surface and volume diffusion processes in crystal growth. J. Cryst. Growth 1971, 8, 79. (42) Ghez, R.; Gilmer, G. H. Analysis of combined surface and colume diffusion processes in crystal-growth: II. Asymmetric capture at steps. J. Cryst. Growth 1974, 21, 93−109.

number of growth units per growth sites

REFERENCES

(1) Vetter, T.; Iggland, M.; Ochsenbein, D. R.; Hänseler, F. S.; Mazzotti, M. Modeling Nucleation, Growth, and Ostwald Ripening in Crystallization Processes: A Comparison between Population Balance and Kinetic Rate Equation. Cryst. Growth Des. 2013, 13, 4890−4905. (2) Xie, B.; Halter, T. J.; Borah, B. M.; Nancollas, G. H. Tracking Amorphous Precursor Formation and Transformation during Induction Stages of Nucleation. Cryst. Growth Des. 2014, 14, 1659− 1665. (3) Tsukimura, K.; Suzuki, M.; Suzuki, Y.; Murakami, T. Kinetic Theory of Crystallization of Nanoparticles. Cryst. Growth Des. 2010, 10, 3596−3607. (4) Wiencek, J. M. New strategies for protein crystal growth. Annu. Rev. Biomed. Eng. 1999, 1, 505−534. (5) McPherson, A. A brief history of protein crystal growth. J. Cryst. Growth 1991, 110, 1−10. (6) Durbin, S. D.; Feher, G. Protein crystallization. Annu. Rev. Phys. Chem. 1996, 47, 171−204. (7) Sleutel, M.; Maes, D.; Wyns, L.; Willaert, R. Kinetic Roughening of Glucose Isomerase Crystals. Cryst. Growth Des. 2008, 8, 4409−4414. (8) Teng, H. H.; Dove, P. M.; De Yoreo, J. J. Kinetics of calcite growth: surface processes and relationships to macroscopic rate laws. Geochim. Cosmochim. Acta 2000, 64, 2255−2266. (9) Nielsen, A. E. Electrolyte crystal-growth mechanisms. J. Cryst. Growth 1984, 67, 289−310. (10) Nielsen, A. E. Mechanisms and rate laws in electrolyte crystalgrowth from aqueous-solution. ACS Symp. Ser. 1987, 323, 600−614. (11) Hellevang, H.; Miri, R.; Haile, B. G. New Insights into the Mechanisms Controlling the Rate of Crystal Growth. Cryst. Growth Des. 2014, 14, 6451−6458. (12) van der Weijden, C. H.; van der Weijden, R. D. Calcite growth: Rate dependence on saturation, on ratios of dissolved calcium and (bi)carbonate and on their complexes. J. Cryst. Growth 2014, 394, 137−144. (13) Wolthers, M.; Nehrke, G.; Gustafsson, J. P.; Van Cappellen, P. Calcite growth kinetics: Modeling the effect of solution stoichiometry. Geochim. Cosmochim. Acta 2012, 77, 121−134. (14) Larsen, K.; Bechgaard, K.; Stipp, S. L. S. The effect of the Ca2+ to activity ratio on spiral growth at the calcite surface. Geochim. Cosmochim. Acta 2010, 74, 2099−2109. (15) Zhang, J.; Nancollas, G. H. Kink Density and Rate of Step Movement during Growth and Dissolution of anABCrystal in a Nonstoichiometric Solution. J. Colloid Interface Sci. 1998, 200, 131− 145. (16) Davies, C. W.; Jones, A. L. The precipitation of silver chloride from aqueous solutions. Part I. Discuss. Faraday Soc. 1949, 5, 103−111. (17) Bracco, J. N.; Grantham, M. C.; Stack, A. G. Calcite Growth Rates As a Function of Aqueous Calcium-to-Carbonate Ratio, Saturation Index, and Inhibitor Concentration: Insight into the Mechanism of Reaction and Poisoning by Strontium. Cryst. Growth Des. 2012, 12, 3540−3548. (18) Stack, A. G.; Grantham, M. C. Growth Rate of Calcite Steps As a Function of Aqueous Calcium-to-Carbonate Ratio: Independent Attachment and Detachment of Calcium and Carbonate Ions. Cryst. Growth Des. 2010, 10, 1409−1413. (19) Burton, W. K.; Cabrera, N.; Frank, F. C. The Growth of Crystals and the Equilibrium Structure of their Surfaces. Philos. Trans. R. Soc., A 1951, 243, 299−358. (20) Bennema, P. Spiral growth and surface roughening: Developments since Burton, Cabrera and Frank. J. Cryst. Growth 1984, 69, 182−197. K



Article

(43) Van Der Eerden, J. P. The advance velocity of steps under the influence of volume and surface diffusion, by direct and indirect incorporation of growth units. J. Cryst. Growth 1982, 56, 174−188. (44) Chernov, A. A. How does the flow within the boundary-layer influence morphological stability of a vicinal face? J. Cryst. Growth 1992, 118, 333−347. (45) Chernov, A. A.; Coriell, S. R.; Murray, B. T. Morphological stability of a vicinal face induced by step flow. J. Cryst. Growth 1993, 132, 405−413. (46) Coriell, S. R.; Chernov, A. A.; Murray, B. T.; McFadden, G. B. Step bunching: generalized kinetics. J. Cryst. Growth 1998, 183, 669− 682. (47) Kandel, D.; Weeks, J. D. Simultaneous bunching and debunching of surface steps - theory and relation to experiments. Phys. Rev. Lett. 1995, 74, 3632−3635. (48) Kandel, D.; Weeks, J. D. Step motion, patterns, and kinetic instabilities on crystal-surfaces. Phys. Rev. Lett. 1994, 72, 1678−1681. (49) Kerisit, S.; Parker, S. C. Free Energy of Adsorption of Water and Metal Ions on the {101̅4} Calcite Surface. J. Am. Chem. Soc. 2004, 126, 10152−10161. (50) Smith, K. S. Strategies to predict metal mobility in surficial mining environments. Rev. Eng. Geol. 2007, 17, 25−45. (51) Parkhurst, D. L.; Appelo, C. A. J. Description of input and examples for PHREEQC Version 3a computer program for speciation, batch-reaction, one-dimensional transport, and inverse geochemical calculations, 2013; p 497. (52) Kowacz, M.; Putnis, C. V.; Putnis, A. The effect of cation:anion ratio in solution on the mechanism of Barite growth at constant supersaturation: Role of the desolvation process on the growth kinetics. Geochim. Cosmochim. Acta 2007, 71, 5168−5179.

L


A Statistical Approach To Explain the Solution Stoichiometry Effect on

Recommend Documents