Phase-Field Modeling of Biomineralization in Mollusks and Corals: Microstructure vs Formation Mechanism

While biological crystallization processes have been studied on the microscale extensively, there is a general lack of models addressing the mesoscale aspects of such phenomena. In this work, we investigate whether the phase-field theory developed in materials’ science for describing complex polycrystalline structures on the mesoscale can be meaningfully adapted to model crystallization in biological systems. We demonstrate the abilities of the phase-field technique by modeling a range of microstructures observed in mollusk shells and coral skeletons, including granular, prismatic, sheet/columnar nacre, and sprinkled spherulitic structures. We also compare two possible micromechanisms of calcification: the classical route, via ion-by-ion addition from a fluid state, and a nonclassical route, crystallization of an amorphous precursor deposited at the solidification front. We show that with an appropriate choice of the model parameters, microstructures similar to those found in biomineralized systems can be obtained along both routes, though the time-scale of the nonclassical route appears to be more realistic. The resemblance of the simulated and natural biominerals suggests that, underneath the immense biological complexity observed in living organisms, the underlying design principles for biological structures may be understood with simple math and simulated by phase-field theory.


INTRODUCTION
Crystalline materials formed by solidification from the liquid state play an essential role in our civilization. 1,2 This class of matter incorporates most of the technical alloys, polymers, minerals, drugs, food products, and so on. Owing to their importance, mathematical models describing the process of crystallization under the respective conditions were and are being developed. Relying on the statistical physical description of phase transitions, the evolving numerical methods, and the ever-increasing computational power, computational materials science reached the level where knowledge-based design of crystalline matter is possible for certain classes of materials (see, e.g., refs 2−4). The models that address the behavior of matter during crystalline solidification range from the molecular time and length scales to the engineering scales. They include ab initio computations; particle-based methods like molecular dynamics (MD), Monte Carlo, or population dynamics simulations and different types of continuum models ranging from the density functional theory of classical particles, via coarse-grained models (such as the time-dependent Ginzburg−Landau, Cahn−Hilliard, and phase-field type order parameter theories that belong to the family of classical field theoretical models widely used in modeling phase transitions of various complexity), to the macroscopic continuum models applicable on engineering timeand length-scales. While this inventory allows the modeling of a substantial range of crystallization phenomena, there are complex cases, for which its use is not straightforward. Such examples are the biomorphic (inorganic) materials 5−10 that form worm-shape or arboresque morphologies by aggregation of crystalline particles, and the process of biomineralization; 11−19 that is, the formation of hierarchically structured organic− inorganic composites in biological systems. Examples of biomineralization include the formation of mollusk shells, 13,14 skeletons of corals 15 and cell walls of diatoms, 16 kidney stones, 17 bones and teeth, 18 and magnetite crystals in the magnetosomes of magnetotactic bacteria, 19 to name a few. The materials formed by biomineralization often have surprisingly good mechanical properties owing to their hierarchical microstructure (see e.g., refs 13,14). Recent imaging and analytic methods provide detailed information on the respective microstructures, which in turn may give clues to the formation mechanism: many of these microstructures are well-known from materials science (such as dendrites, spherulites, cellular, columnar shapes, etc.). 11−15,17−19 This raises the possibility that with some adjustment/further development, the models developed in materials science can be used to reverse engineer the biomineralization process and learn the pathways used by nature to create these complex structures, which may inspire new technologies for creating novel composite materials. 20−25 Recently, we explored the possibility of developing predictive mathematical models for biomorphic crystallization and for relatively simple biomineralization processes by adopting wellestablished methods of computational materials science and adjusting them to the circumstances as necessary. 26 −28 The research done so far is confined yet to relatively simple cases of extracellular biomineralization such as mollusk shell formation 26,27 or microstructure evolution of spherulitic structures in coral skeletons 28 but is expected to deepen the general understanding in the field, and the tools developed in the course of this research might open the way for modeling more complex cases of crystallization in biological systems such as formation of bones, kidney stones, and so on.
In the present paper, we concentrate on the modeling aspects of such an approach, outlining possible minimum requirements for phase-field modeling of biological crystallization processes, and demonstrate that with appropriate choice of the model parameters and boundary conditions phase-field models can approximate the polycrystalline microstructure formed in simple cases of biomineralization (shell formation in mollusks such as Figure 1. Schematic view of bivalve molluscan anatomy with successive magnification of the mantle−nacre interface. 39 A thin liquid-filled extrapallial space is indicated (its thickness and content is open to debate). The interlamellar membrane is made of a viscoelastic chitin-based organic substance, whereas the mineral constituent is crystalline CC (aragonite). (Reproduced with permission from ref 39. Copyright 2009 United States National Academy of Sciences). Here "growth direction" indicates the lateral growth of the shell. However, this is based on layerwise growth perpendicular to the surface via formation of new layers that also grow laterally, a process which is responsible for both the thickening and the sidewise spreading of the shell. In the present work, we model the local thickening of the shell. bivalves, gastropods and cephalopods, and sprinkle formation in coral skeletons).

Microstructures Formed during Biomineralization
Before outlining the phase-field models we used in the present research, we give a short account of the experimental results on the observed microstructures. Polycrystalline microstructures formed in biomineralization processes have been investigated by a variety of experimental methods, including optical microscopy (OM), scanning electron microscopy (SEM), electron back scattering diffraction (EBSD), X-ray tomography, and polarization-dependent imaging contrast (PIC, 29,30 ), among others. Here we concentrate on experimental results obtained on the microstructure of molluscan shells and coral skeletons.
1.1.2. Mollusk Shells. Mollusk shells are complex organomineral biocomposites with a broad range of species-dependent microstructures. 13,14,30−39 A schematic view of bivalve anatomy having a nacre-prismatic shell is shown in Figure 1. Moving inward, the sequence of the individual layers is as follows: the organic periostracum, a leathery "skin", that encloses the domain where biominerlization takes place, a dominantly mineral (calcium carbonate, CC) prismatic layer, and the nacre composed of CC tablets and organic interlamellar membranes, the submicron thick extrapallial liquid, 31,34−36 and then the outer calcifying epithelum layer of the mantle. Images showing typical microstructures of mollusk shells of similar type are displayed in Figure 2. 13,26,37,38 A recent study shows that in members of three classes of mollusks Unio pictorum (bivalve), Nautilus pompilius (cephalopod), and Haliotis asinine (gastropod), the shell displays a common sequence of ultrastructures: a granular domain composed of randomly oriented crystallites, a prismatic domain of columnar crystallites, and the nacre 26,27 (Figure 3). It has been shown that the layered structure of nacre may contain screw dislocation-like defects (see Figure 4). 14,39−44 Of these structures, the prismatic layers show mechanical flexibility, whereas the nacre (also called "mother of pearl") is fairly rigid but hard; the combination of the two yields a surprisingly strong yet flexible biocomposite (see e.g., ref 22  45 and the connecting living tissue, which secretes calcium carbonate to create a hard shelter (the corallite, a tubular hollow structure on which the polyp sits), 46 into which the polyp can retreat if danger is detected (see Figure 5a,b). The polyps are transparent, their color originates from photosynthesizing algae (zooxanthallea) that live in symbiosis with the polyp and feed the polyp sugars and oxygen. The surface of the skeleton is intricately structured, 46 depressions, ridges, cavities are arranged into complex patterns reflecting the radial symmetry of the polyp (Figure 5b). The CC crystal (aragonite) building the porous skeleton ( Figure 5c) has a spherulitic microstructure, that is, a radial arrangement of crystals radiating from a common center. The center does not have to be a point; it can be a line or a plane in plumose spherulites. In the case of coral skeletons, the centers are curved planes, termed centers of calcification (CoCs) (Figure 5d). 15,47 From these CoCs, acicular fibers grow radially and then arrange into fan-like bundles that finally group into a feather-duster-like shape termed "trabecula". 15 Besides the plumose spherulitic structure, randomly oriented nanoscale crystallites "sprinkles" are also present 28,48,49 (Figure 5e). Recent PIC mapping experiments performed at synchrotron indicate that the amount of the  submicron-sized sprinkle crystallites varies considerably among different coral species: in some of them, the sprinkles are missing (Figure 5f), whereas in others, submicron size crystallites appear at the perimeter of the skeleton, along grain boundaries, at the growth front of spherulitic trabeculae, or in bands formed in the interior of the skeleton. 28 The discovery of sprinkles inspired a refined model for spherulitic growth in corals, described in detail in ref 28. Randomly oriented sprinkles are the first nucleated crystals at the growth front. With further growth, those oriented radially have space and thus continue to grow, and those oriented tangentially run into each other and stop growing. This is why they stay small. A coarsening process then makes the larger crystal grow larger at the expense of the smaller ones, which disappear. In most mature spherulites, therefore, no sprinkles remain. In the skeleton of some coral species, however, some sprinkles do not disappear, presumably because they are kinetically trapped.
In our previous work, we modeled this process by phase-field simulations 28 and raised the possibility that other spherulites may grow this way, including aspirin, chocolate, and geologic crystals. However, the formation mechanism of sprinkle bands and the origin of different amount of sprinkles in the skeleton of different coral species is not yet fully understood. Mathematical modeling is expected to help to identify the governing factors. Molecular/ionic mobility in the calcifying fluid is expected to be orders of magnitude higher than in the solid. As a result, in the case of growth via ion-by-ion addition, either a slow supply of the ions or a kinetic barrier of ion deposition can keep the growth rate sufficiently low. This offers limitations to the possible mechanisms, as will be discussed later.
1.1.4. Biomineralization on the Nano-and Macroscale. The complexity of biomineralization stems mainly from the fact that the fluid and solids incorporate organic molecules, the role of which is largely unknown. 35,36,50,51 For example, nanoscale amorphous calcium carbonate (ACC) globules are essential for the formation of mollusk shells and coral skeletons. 13,14,52−54 Possible pathways from free ions to crystalline CC are reviewed in ref 55. Evidently, all the possible processes cannot be explicitly incorporated into an orientation field based phase-field (OF-PF) approach designed to model polycrystalline microstructures on the mesoscale. In this work, we investigate two specific cases: (i) formation of crystalline CC via the classical route of ion-by-ion addition of Ca 2+ and CO 3 2− , and (ii) the nonclassical route via an amorphous precursor (ACC) deposited at the solidification front. In the latter case, we hypothesize that the crystallization rate is determined by the velocity of crystal growth into the ACC layer and not by the rate of supplementing ACC; that is, growth is controlled by the selfdiffusion in ACC. This hypothesis can account for the typically months' time scale of shell/skeleton formation, which would be difficult to interpret, for example, on the basis of ion-by-ion deposition directly from the extrapallial fluid or other aqueous solutions. A specific realization of mechanism (ii) is presented in ref 56. The extrapallial fluid contains various ions and organic molecules as shown by in vivo studies. 35, 36 1.1.5. Growth Rate of Mollusk Shells and Coral Skeletons. The shell of bivalves grows typically by 100−300 μm lunar-day increments, 57,58 corresponding to a growth rate of about v ≈ 4.2 × 10 −11 −1.3 × 10 −10 m/s, which decreases with age. 57 The thickening rate of the shell of Tridacna deresa was estimated to be v ≈ 1.6 × 10 −10 −4.9 × 10 −10 m/s in its early life, which decreases to v ≈ 3.2 × 10 −11 −2.3 × 10 −10 m/s in the later life. 59 Comparable growth rates were reported for freshwater gastropods v ≈ 9.5 × 10 −11 m/s. 60 In contrast, the coral skeleton growth rates range between about 1 and 37 cm/year, corresponding to 3.2 × 10 −10 to 1.2 × 10 −8 m/s. 61−63

MODELING SECTION
There are two main categories of the phase-field (PF) models developed to address polycrystalline freezing: (a) the multiphase-field (MPF) models that assign a separate phase field for every crystal grain 64−70 and (b) the orientation-field based phase-field (OF-PF) approaches, in which the local crystallographic orientation is monitored by a scalar field (2D) 71−78 or quaternion/rotation matrix fields (3D). 78−85 Recent developments in these areas were reviewed in refs 70 and 78, respectively.
Both the MPF and OF-PF approaches have their advantages, yet complex polycrystalline growth forms (such as disordered dendrites, crystal sheaves, various types of spherulites, and fractal-like polycrystalline aggregates, etc.) were so far modeled exclusively using the OF-PF approach. A further advantage of these models is that they allow a continuous variation of the orientation, a feature particularly useful in modeling biomineralization.
In the OF-PF models, the local phase state of the matter is characterized by a coarse grained structural order parameter, the "phase field", that is time-and space-dependent and monitors the transition between the liquid and solid states. This field is usually coupled to other slowly changing fields, such as the concentration field of the constituents, and the orientation field. The free energy contains bulk free energy and gradient terms for these fields. The equations of motion can be obtained via a variational route. The thermal fluctuations are represented by adding noise (nonconserved fields) or flux noise (conserved fields) to the equations of motion that satisfy the fluctuation− dissipation theorem, allowing for homogeneous nucleation 72 The inclusion of foreign surfaces/particles, represented by appropriate boundary conditions that set the wetting properties, allows the modeling of heterogeneous nucleation 83−85 and solidification in confined space. 77,86 The addition of phase-field noise makes the modeling of Brownian motion of solid particles possible. 87 The OF-PF models may also be coupled to hydrodynamic flow, which can be represented by either the Navier−Stokes equation 88,89 or the lattice Boltzmann technique. With an appropriate overlapping grid technique, flow coupled motion of growing solid particles can also be treated. 90 The main virtue of PF modeling is that the growth morphology can be computed on the basis of the freezing conditions and the thermophysical properties. Images obtained by OF-PF modeling in two-and three dimensions (2D and 3D) 74,78,84,86 are compared with experiments 74,86,91−101 in Figure 6.
Recently, PF modeling has been extended for microstructure formation in mollusk shells and coral skeletons, and promising agreement was seen between experiment and the predicted microstructures. 26−28 We give below a detailed account of the modeling efforts and present new results.
Herein, we continue further the quest for a minimum phasef ield model of the biomineralization process during the formation of mollusk shells and coral skeletons. Evidently, we cannot model the living organism in the framework of this approach; their functions are represented by appropriate boundary conditions. Unfortunately, usually little is known of the thermodynamics of the multicomponent fluids involved, of the interfacial free energies and their anisotropies, and the respective diffusion coefficients. Therefore, we aim at identifying the main components that a minimal theory needs to incorporate for qualitatively reproducing the microstructures seen in the experiments.
In the present study, we rely on three specific formulations of the phase-field theory. In two of them, 26,27 the local state is characterized by the phase field, ϕ(r, t) that monitors the process of crystallization, the concentration field c(r, t) specifying the local composition, and a scalar orientation field θ (r, t), which represents the local crystallographic orientation in 2D. The latter field is made to fluctuate in time and space in the liquid, a feature that represents the short-range order present in the liquid state (as done in refs 72,74−79 and 83−86). The third model does not include an orientation field. It has been developed to handle two-phase spiraling dendrites during eutectic solidification in ternary alloys in 3D, 102 assuming a fixed relative orientation of the solid phases.
In Figure 7, we summarize the problems addressed herein, the models used, the assumed micromechanisms, and the corresponding phase transitions.
The first model will be used to address four cases of biomineralization: (i) mollusk shell formation from aqueous solution by ion-by-ion addition (ions from the extrapallial fluid attach to the surface of the crystalline phase), (ii) mollusk shell formation via amorphous precursor (ACC → CCC transition), (iii) formation of columnar nacre via ion-by-ion attachment, and (iv) formation of coral skeletons via ion-by-ion attachment. In all these cases, the same mathematical model will be used, however, with specific input parameters and initial-and boundary conditions. Model PF1 is defined by eqs 1−6 shown below. This model is similar to the standard binary PF theory by Warren and Boettinger; 103 however, it is supplemented with an orientation field as done in refs 72,74−78 Accordingly, the free energy of the heterogeneous system is expressed as where T is the temperature. Parameters ε ϕ are expressed in terms of the free energy γ i , the thickness δ i of the crystal−liquid interface, and the melting point T i of the i th pure component (i = A or B that stand for the organic and mineral components, respectively). In model PF1 ε c 2 = 0 is chosen. s = s(ϑ, θ) = 1 + s 0 cos{kϑ − 2πθ} is an anisotropy function corresponding to an interfacial free energy of k-fold symmetry and strength s 0 , whereas ϑ = arctan(ϕ y /ϕ x ) is the angle of the normal of the interface in the laboratory frame, while ∇ϕ = [ϕ x , ϕ y ]. The angular (circular) variables ϑ and θ are normalized so that they vary between 0 and 1. The bulk free energy density reads as, and varies between the free energy densities of the crystal and mother phases (f C and f M , respectively) as prescribed by the interpolation function p(ϕ) = ϕ 3 (10−15ϕ + 6ϕ 2 ). Here f C and f M were taken from the ideal solution model. The orientation free energy density is as follows: where the parameter H can be used to tune the magnitude of the grain boundary energy. Owing to the scalar orientation field, model PF1 is applicable exclusively in 2D. The time evolution of the heterogeneous system is described by variational equations of motion (EOMs): where M ϕ , M c , and M θ are mobilities that determine the time scale of the evolution of the individual fields, and are related to coefficients of the self-diffusion, interdiffusion, and rotational diffusion. 72,74−78 The chemical and orientation mobilities are made phase-dependent as where i = c or θ, and indices M and C denote values for the mother and crystalline phases. The corresponding dimensionless mobilities are defined as Here v m is the average molar volume of the components, R the gas constant, ξ the length scale, whereas HT is the energy scale of the grain boundary energy. Gaussian white noise terms ζ i are added to the EOMs to represent the thermal fluctuations (here i = ϕ, c, and θ). 3D generalizations of model PF1 can be found elsewhere, 78−82 which, however, require quaternion or rotation matrix representation of the crystallographic orientation, as opposed to the scalar field used here.

Phase-Field Model 2 (PF2)
The second model will be used to provide a refined model of mollusk shell formation including the nacreous structures in the case of (i) ion-by-ion attachment and (ii) amorphous precursor mediated process. This OF-PF model was originally developed to describe eutectic solidification, while keeping a fixed orientational relationship between the two solid phases inside the crystal grains. 104 To realize this, the square-gradient term, (1/2)ε c 2 T (∇c) 2 , was retained in the free energy density (choosing ε c 2 = 2ε ϕ 2 ), and a more complex orientational free energy term was used that realizes a fixed orientational relationship at the solid−solid phase boundaries. Here h(  . Summary of biomineralization-related problems that are addressed within this work using phase-field methods. In each case, the applied model, the respective micromechanism, and the relevant phase transition are specified. Phase-field models PF1, PF2, and PF3 are defined in the text.

Phase-Field Model 3 (PF3)
This model will be used to address the formation of 3D topological defects in sheet nacre via crystal growth into hACC (here the third component plays the role of water that has smaller solubility in CCC than hACC). Modeling of sheet nacre by PF3 is a 3D analogue of the amorphous precursor mediated case in PF2. The free energy of the crystallizing system is given by the expression 102 where c = [c 1 ,c 2 ,c 3 ] and Σ i c i = 1, whereas w and ε c 2 are constants. The bulk free energies of the solid and liquid phases are taken from the regular and ideal solution model: and where Ω i,j are the binary interaction coefficients in the solid. The respective EOMs obtained variationally are as follows and where M c is the 3 × 3 mobility matrix. With the choice of 1 and −0.5 for the diagonal and off-diagonal elements, the criterion Σ i c i = 1 is automatically satisfied. The diffusion is switched off in the bulk solid. Further information on model PF3 is given in ref 102.

Numerical Solutions
The EOMs were solved numerically in a dimensionless form, on rectangular uniform grids, using finite difference discretization with forward Euler time stepping. The PF1 and PF2 codes were run on a CPU cluster of 480 CPU cores using MPI protocol. Typical runs on a 1000 × 2000 grid took between about 8 to 15 h, depending on the number of time steps that varied from 2.4 × 10 5 to 5 × 10 5 , as required by the velocity of crystallization. The code for PF3 was run on high-end graphics processing units (GPUs) and was solved on a 3D rectangular grid.

Materials' Parameters
We review here the present status of input data required for a quantitative modeling. The OF-PF models require a fairly detailed information on the systems studied. This incorporates the free energy of all the relevant phases as a function of temperature and composition; all the interface energies; and the translational, chemical, and rotational diffusion coefficients.
Since in the biomineralization problems, we address here, the dominant CC polymorph is the metastable aragonite, we use the thermophysical data available for this polymorph, as much as possible.
In case, where no information is available, we use values for another CC polymorph (calcite). The input data are collected in Tables 1 and 2 for Model PF1, in Tables 3 and 4 for Model PF2, a few common ones are presented in Table 5. In these Tables, "aq. sol. → CCC" indicates data relevant to the ion-by-ion mechanism, whereas "ACC → CCC" denote those for the amorphous precursor mediated case. 2.5.1. Thermodynamics. Unfortunately, only limited thermophysical information is available even for the pure CC system from experiment and MD simulations, such as the phase diagrams 105,106 and the equilibrium shapes reflecting the anisotropy of the interface energy. 107,108 During biomineralization, however, a variety of ions and organic macromolecules are present that may influence/control the crystallization process. 109−115 Accordingly, it is a nontrivial task to obtain accurate input data for mesoscale modeling; for example, selective adsorption of ions or organic molecules on different crystal faces may change growth morphology 111,112 or influence the formation of polymorphs of CC. 114,115 Owing to this lack of information, we present here generic approaches that are based on simplified hypothetical model systems of properties similar to those used in refs 26−28.
2.5.2. Diffusion Coefficients. As noted above, we address here two scenarios for the formation of crystalline CC (CCC): (i) diffusion controlled growth of crystalline CC directly from The subscripts M and C stand for the mother and crystalline phases. The chemical mobility of the former was used as reference, as its chemical diffusion coefficient was used in making the EOMs dimensionless.   The subscripts M and C stand for the mother and crystalline phases. The chemical mobility of the former was used as reference, as its chemical diffusion coefficient was used in making the EOMs dimensionless.

JACS Au
pubs.acs.org/jacsau Article aqueous solution via ion-by-ion addition; (ii) diffusion controlled growth of crystalline CC into a hydrated ACC (hACC) layer that is assumed to be deposited on the solidification front (by vesicles or the-ion-by-ion process) with a sufficient rate so that deposition is not the rate-limiting process. These scenarios differ in the diffusion coefficient we assign to the "mother phase" (either aqueous solution or ACC) that crystallizes. Since the equations of motion were made dimensionless using the diffusion coefficient of the mother phase, and we assume that the relative magnitudes of M ϕ , M c , and M θ remain the same in the mother phase independently whether it is liquid or amorphous, the two scenarios differ in only the dimensionless mobilities assigned to the crystalline phase. Aqueous Solutions. The coefficient of ion diffusion in aqueous solutions at room temperature is in the order of D L ≈10 −9 m 2 /s. 116 Amorphous CC. The ion diffusion in ACC at 300 K is in the order of D ACC,ion ≈ 10 −15 m 2 /s. 117 However, the diffusion coefficient of the water molecules from MD simulations is typically D hACC,H 2 O ≈ 10 −14 −10 −13 m 2 /s for the slow H 2 O molecules, although a few percent of H 2 O molecules that have orders of magnitude faster diffusion (D hACC,H 2 O ≈ 10 −11 m 2 /s) are also present. 117,118 However, biogenic ACC is almost anhydrous. 119 Therefore, water diffusion is expected to play a negligible role. The rate limiting factor for the structural transition is expected to be the slowest of these processes; accordingly, we use the diffusion coefficient for the ions, D ACC,ion ≈ 10 −15 m 2 /s. 117 Crystalline CC. We are unaware of self-diffusion data for aragonite. There are, however, experimental data for calcite. The Mg diffusion coefficient in calcite is about D calcite ≈ 10 −21 m 2 /s at 823 K, which, extrapolates to D calcite ≈ 10 −53 m 2 /s at room temperature provided that the diffusion mechanism does not change. 120 A different estimate is obtained via extrapolating the Mg−Ca interdiffusion data of ref 121 to room temperature, which yields D calcite ≈ 10 −29 m 2 /s at 300 K for Mg diffusion in calcite. An even higher value emerges from the radioactive tracer method D calcite ≈ 10 −23 m 2 /s. 122 Summarizing, herein, we opt for D ACC,ion ≈ 10 −15 m 2 /s and D calcite ≈ 10 −29 m 2 /s for ion diffusion in ACC and calcite, respectively; assuming thus that the diffusion data for calcite can be viewed as a reasonable order-of-magnitude estimate for the other polymorphs, including aragonite. The corresponding dimensionless mobility data are presented in Tables 1 and 3. The data of magnitude 10 −20 to 10 −12 indicate that in the crystal, independently of the mother phase, the reduced mobilities are essentially zero for the concentration and orientation fields; that is, these fields do not show time evolution within the crystal on the time scale of the simulations.
2.5.3. Interfacial Free Energies. The experimental and theoretical results for the water-aragonite interfacial free energy are about 150 mJ/m 2 . 123,124 A recent ab initio theoretical treatment provides a considerably larger value (280 mJ/m 2 ), and information on its anisotropy for small aragonite clusters. 125 Herein, we use 150 mJ/m 2 . We are unaware of data for the free energy of the ACC-aragonite interface. Using Turnbull's relationship 126 for the interfacial free energy, γ = αΔH/ (N 0 v mc 2 ) 1/3 , where α is a constant, ΔH heat of transformation, N 0 the Avogadro-number, and v mc is the molar volume of the crystalline phase, a crude estimate can be made on the basis of the enthalpy difference between ACC and calcite: ΔH calcite-ACC = (14.3 ± 1.0) kJ/mol. 119 A similar value may be expected for aragonite. 119 Considering α = 0.55 from MD simulations, 127 one obtains γ aragonite-ACC ≈ 87 mJ/m 2 , which result, however, needs independent confirmation by other experimental/theoretical methods.
Once the thermodynamic data are fixed for components A and B, and the interfacial free energy is given for one of the components, models PF1 and PF2 predict the interfacial free energy for the other component, provided that the interface thicknesses are similar. 103 For materials of comparable entropy of transformation, this realizes γ ∝ T trans , where T trans is the temperature of the phase transition, a relationship that works well for the solid−liquid interfacial free energy of metals. 128 2.5.4. Qualitative Modeling. Despite our efforts to collect a full set of the required materials parameters, owing to uncertainties of the thermodynamic diving force of crystallization and of the interface energy estimates, the simulations we present can only be regarded as qualitative. They are aimed at demonstrating that phase-field modeling has the potential to capture various microstructural/morphological aspects of biomineralization. This summary of the present status of input data may give hints where further experiments and microscopic theory can help mesoscale modeling.

Modeling of Microstructures Mimicking Mollusk Shells
Herein, we hypothesize that (1) the granular domain is produced by heterogeneous nucleation (either on the surface of periostracum or on organic particles); and (2) the physical background of the columnar → nacre morphological transition is the observation that decreasing the driving force of solidification, an initially diffusionless process (full solute trapping) is replaced by partitioning, which first appears in the form of alternating layers rich in one or the other component (see Figure 8). 129,130 Here Ω M,C = Ω 0,M,C − Ω 1,M,C T. The latter mechanism is the one that causes band formation in Figures 9 and 10. Further assumptions made during the application of these models (e.g., boundary conditions specific to the studied problems) are recapitulated below. In the case of PF1, reduction of the driving force leads to a transition of a chemically homogeneous solid to alternating solid−liquid layers, whereas in the case of PF2 partitioning appears via alternating α−β bands. For even smaller driving forces, PF1 would produce seaweed/dendritic structures, whereas PF2 would yield lamellae perpendicular to the growth front.
3.1.1. Shell-Like Microstructure in Model PF1. In a recent OF-PF study, 26 we made the following assumptions that define the conditions under which eqs 1−6 were solved, when modeling the formation of mollusk shells within model PF1: • The CC crystals grow into the extrapallial fluid by the molecule/ion attachment mechanism. • Binary ideal solution thermodynamics (CC and organic component) is applied. Evidently, treating the extrapallial fluid as a quasi-binary solution is a gross simplification. During crystallization of the CC-rich crystal, and an organic-component-rich "fluid" forms from the original homogeneous mixture. This construction was used as a simple means to provide thermodynamic driving force for CCC precipitation. In the present work, besides this, we explore a different scenario shown in Figure 7, in which the CCC crystal grows into an ACC precursor that forms continuously ahead of the crystallization front. The CCC front propagates into this ACC

JACS Au
pubs.acs.org/jacsau Article layer and has no direct contact with the extrapallial fluid. It is also assumed that the formation of the ACC layer is fast enough, so that it is not the rate-limiting process. The respective amorphous → crystal transition has in principle a reduced driving force, while orders of magnitude smaller diffusion coefficients prevail in the amorphous phase, yielding a far longer time scale for crystallization, when compared with crystallization from an aqueous solution. In this scenario, the phase field monitors the amorphous → crystal transition, rather than the crystallization of a liquid. We furthermore assume that the coefficients of the translational, chemical, and rotational diffusion decrease proportionally by orders of magnitude during this process, retaining the same relative magnitudes of the mobilities as in the liquid (see Table 1).
Since the EOMs are solved in dimensionless form, where time is dedimensionalized using the chemical diffusion coefficient of the mother phase, the dimensionless chemical mobilities of the mother phase remain unchanged. What differs between the present computation and the previous one in ref 26 is the magnitude of the individual mobilities in the crystalline and the mother phase (see Table 1). Following the general principles of phase-field modeling (see e.g. ref 131), the phase-field mobility (b, f) columnar growth via directional solidification in concentration gradient yielding the prismatic structure, and (c, g) layerwise formation of alternating CCC and organic layers (sheet nacre). In panels (a−c) the mother phase is an aqueous solution (D c,M = 10 −9 m 2 /s) in the simulation shown; whereas an amorphous precursor (D c,M = 10 −15 m 2 /s) is assumed in the simulation shown in panels (d−f). The wavelength of the layered structure is roughly proportional to the free energy of the mother phase − CCC interface. The thickness of the mother phase (extrapallial domain) is assumed to be constant. In (a−c) and (e−g) different colors stand for different crystallographic orientation, while colors white and black stand for the mantle of the mollusk and the mother phase, respectively. In (d) and (h), gray, yellow, and blue indicate the mother phase, the organic phase, and CCC. The qualitative phase-field simulations were performed on 2000 × 200 grids (corresponding to an area of 26.25 μm × 2.625 μm with the present choice of model parameters). In (a−c) and (e−g) time elapses downward, whereas (c) and (d) and (g) and (h) display snapshots taken at the same time.

JACS Au
pubs.acs.org/jacsau Article is assumed to be independent of the phase. In contrast, in agreement with the experimental diffusion coefficients, the chemical and orientational mobilities are assumed to be 20 orders of magnitude larger in the extrapallial fluid then in the crystal, and 14 orders of magnitude larger in ACC than in CCC. We retain the assumptions made in ref 26, as listed above with the difference that in the nonclassical mechanism it is the CC content of the ACC layer that decreases exponentially (while the organic content increases) with the distance from the periostracum. For the sake of simplicity, we employ the same dimensionless driving force as in ref 26 for both the classical and the nonclassical cases.
The microstructures that evolved in the two cases are compared in Figure 9. The characteristic microstructural transitions are present in both simulations. Whether from the fluid or the amorphous phase, first small randomly oriented CCC grains form in the neighborhood of the periostracum, of which crystal grains grow further inward yielding elongated crystals of random orientation that compete with each other. With increasing distance from the periostracum, that is, with decreasing supersaturation, growth slows down, and the separation of the two constituents becomes possible, which leads to the formation of alternating CCC and organic-rich layers, closely resembling the nacre. The mineral layers in the "nacre" are composed of segments of different orientations (see Figure 9b,d), which can be viewed as a 2D analogue of the usual 3D mineral platelets. The predicted sequence of the morphological transitions is similar for both mechanisms (i.e., for aq. sol. → CCC and ACC → CCC); however, there are minor differences in the relative thicknesses of the granular, columnar prismatic, and layered nacre structures.
While the respective microstructures are rather similar, the typical size scales for the ACC → CCC transition is smaller than for the aq. sol. → CCC, roughly proportionally with the respective interfacial free energies. Alternating CCC and organic-rich layers akin to sheet nacre form here due to a process described in ref 129, with the difference that under the present conditions a roughly flat "banded structure" forms, and thermal diffusion is replaced by chemical diffusion. Apparently, the orientational information is only partly transferred through the organic layers. Mineral bridges (discontinuities of the organic layers) are also observed. The granular → columnar → layered morphological transitions occur here because of changes in the growth velocity.
At high supersaturations nucleation and dif f usionless solidification takes place forming the granular domain. At medium supersaturations nucleation ceases, only competing growth of the existing particles takes place yielding the columnar domain, whereas at small supersaturations alternating CCC and organic layers occur, forming a layered structure that closely resembles the sheet nacre.
Diffusionless crystallization is possible, when the diffusion length l D = D c,M /v is comparable to the thickness of liquid/ crystal or amorphous/crystal interfaces (d ≈ 10 −9 −10 −8 m 132,133 ), where v is the growth rate. The transition from diffusion controlled to diffusionless growth takes place in the regimes 10 −4 < dv/D c,M < 1 or 10 −2 < dv/D c,M < 10, depending on the model. 134 Considering l D ≈ 10 −9 m and a typical experimental growth rate of v ≈ 10 −10 m/s, one finds that the mother phase needs to have a diffusion coefficient of D c,M ≈ 10 −15 −10 −20 m 2 /s to show transition toward diffusionless growth on the time scale required. This clearly rules out the possibility that the CC crystals form dominantly by direct ion-by-ion attachment from the extrapallial fluid as D c,M ≈ 10 −9 m 2 /s applies for the latter process. In turn, this range of D c,M is consistent with crystallization from amorphous CC, as the magnitude of D c,M falls in the range diffusion coefficient takes in the amorphous state. 135 Summarizing, while the two mechanisms considered here lead to similar microstructures, crystallization via the ACC precursor seems preferable to direct solidification via ion-by-ion addition from the extrapallial fluid, as in the latter case the diffusion coefficients are rather high, that is, crystallization is expected to be fast, unless the "fluid" is of high viscosity (not realistic for the extrapallial fluid). Another problem of direct precipitation from the aqueous solution is that due to the high diffusion coefficient, the assumed initial exponential spatial dependence of CC supersaturation is only temporary on the time scale of shell growth, unless crystal growth is so fast that diffusional equilibration cannot take place. This is, however, at odds with the experimental growth rates.
Note that modeling of the experimentally observed orientational ordering in the columnar (prismatic) and nacre domains that yields coalignment of the c' axis of the crystallites, requires a 3D orientation field. Work is underway in this direction.
Finally, we note that the thickness of the organic and mineral layers in the nacre are typically 10−40 nm and 300−600 nm, 136−140 respectively. Our qualitative simulations give a considerably larger relative thickness for the organic layer. This is partly because we intended to model the whole granular → columnar → nacre sequence, and because of limitations of available computational power, we cannot have sufficient spatial resolution to realize a more realistic thickness ratio. If modeling is limited to the nacre, one is expected to achieve a better agreement.
3.1.2. Shell-Like Microstructures in Model PF2. To relax some of the simplifying assumptions made in PF1, a refined model (PF2) was proposed for modeling the formation of mollusk shells in ref 27. In this model, two solid phases form simultaneously from the liquid state, a mineral-rich and an organic-rich, while a fixed relative orientational relationship is forced between the solid phases formed inside the same crystal grain. This realizes a strong orientational coupling between the solid phases. Simultaneous formation of two solids occurs, for example, in eutectic or peritectic systems. In the refined approach, we opted for the former case. The main assumptions that set the conditions, under which eqs 1, 2, 4, 6, and 7 defining model FP2 were solved, are as follows: 27 • Crystal growth of CC happens via molecule/ion attachment. • A binary eutectic model thermodynamics (regular solution) applies. • The mineral content of the extrapallial fluid emitted at the surface of the mantle decreases exponentially with time. • Formation of granular CC crystals starts by heterogeneous nucleation on organic heterogeneities, whose number density is assumed to decrease exponentially with the distance from the periostracum. • The thickness of the extrapallial domain (distance between the mantle and the solidification front) remains constant. (In the simulation, the position of the mantle surface varies in accord with the solidification rate.) In this approach, the assumption that CC supersaturation decreases toward the mantle is removed and is replaced by the more natural assumption that the CC supersaturation at the JACS Au pubs.acs.org/jacsau Article mantle decreases exponentially with time that leads to an analogous result, however, in a more natural way. The simulation presented in ref 27 shows a good qualitative agreement with the experimental microstructures observed in various types of mollusk shells. It is reassuring that model PF2 recovers the experimental microstructure though in a somewhat more ordered form, which could however be made more random by varying the noise added to the equations of motion. At the present state of affairs, it is difficult to decide which of models PF1 or PF2 should be considered superior to the other, yet the two-solid model is probably closer to the reality. Herein, we explore whether the microstructure remains similar, when assuming ACC mediated crystallization in the framework of model PF2. The results obtained by the two mechanisms of crystallization are compared Figure 10. Apparently, as in the case of model PF1, the two micromechanisms for crystalline CC formation yield similar microstructures.
Note that in model PF2, the nacre is composed of two alternating solid phases. Under the conditions used in our simulation the organic component remains in an amorphous or nanocrystalline state. Remarkably, the predicted "nacre" structure recovers such details of the experimentally observed microstructure as the "mineral bridges" and the line defects across the organic layers termed "aligned holes". 33,34,138,140 We will show, however, in the next section that despite this close similarity, the formation mechanism of the nacre can be considerably more complex than predicted here.
3.1.3. Discussion of Results from Models PF1 and PF2. First, we compare the solidification rates obtained from the simulations performed assuming D L ≈10 −9 m 2 /s (aqueous solution) for the mother phase, with those using D ACC,ion ≈ 10 −15 m 2 /s (ionic diffusion in ACC). We wish to stress that the velocities evaluated from the simulations can only be considered qualitative, as the thermodynamic driving force we used in this study might be well away from the true ones, which in turn may differ from the value obtained for the pure aqueous solution− CCC system, or the ACC−CCC system. Work is underway to perform more quantitative phase-field simulations to determine the growth and dissolution rates in pure systems, for which reasonably accurate experimental data are available.
The growth rate results are presented in Figure 11. As one expects the average growth rate is roughly proportional to the diffusion coefficient of the mother phase. The growth velocities predicted for crystallization from ACC by models PF1 and PF2 are about two and one orders of magnitude higher (v ≈ 10 −8 and 10 −9 m/s, respectively) than the experimental ones 57−60 (≈10 −10 m/s, see gray zone in Figure 11), whereas the rates predicted for the ion-by-ion addition are about v ≈ 10 −1 and 10 −2 m/s, which are about 8−9 orders of magnitude too high. On this ground, the mechanism based on the fast ion-by-ion addition can be excluded. One may perhaps argue that the production of calcium and carbonate ions in the surface layer of the mantle (outer epithelium) may be the rate-limiting process, which may be then taken so slow as to match the experimental growth rate. However, in that case, it is not the diffusion in the mother phase that controls the time scale of the process, and thus, the mechanisms models PF1 and PF2 rely on would not be present.
The velocity vs time relationships show characteristic differences for the two models: PF1predicts a steeply increasing velocity, followed by a plateau decreasing slowly with time for the domain of alternating solid−liquid layers. In the present simulations, the early stage behavior of the two models is different due to reasons different of the growth mode: while model PF1 starts with surface induced heterogeneous nucleation on the periostracum, PF2 relies on volumetric heterogeneous nucleation on organic impurities (note that both heterogeneous nucleation mechanisms can be adopted in both models). As a result, in model PF2, fast initial crystallization is observed during the formation of the granular layer via volumetric heterogeneous nucleation. This is followed by steady-state growth (roughly constant growth velocity) in both the prismatic and the nacreous domains, due to the lack of long-range diffusion during fast eutectic solidification. In contrast, in PF1 a continuously decelerating solidification is observed, which is combined with oscillating growth rate in the nacre. The hypothesized mechanism (volumetric heterogeneous nucleation in a thick layer at the periostracum) that creates protection for the mollusk fast in the early stage of shell formation is advantageous from the viewpoint of survival.
We note that because of computational limitations the maximum linear size of the computational domain we used was about 26.4 μm both in the PF1 and PF2 simulations. However, analogous structures can be produced on a larger size scale via reducing the rate by which the driving force of crystallization decreases with position/time and with an extended initial domain filled with heterogeneous nucleation centers.
It is appropriate to mention that while our models describe the formation of the granular and prismatic layers and the sheet nacre reasonably well within the framework of directional solidification, the predicted mechanism for the formation of the nacre via alternately precipitating mineral and organic layers may be oversimplified. This is especially true for columnar nacre: experiments indicate that the formation of a quasi-periodic network of organic membranes precedes the formation of the CCC layers, which fill the space between the organic membranes later, as illustrated in Figure 12. 136−140 Although one could imagine that the organic membranes form periodically in space via an oscillating chemical reaction of the diffusion-reaction type, while the extrapallial fluid fills the space between the membranes, apparently, the real mechanism is more complex: first a multilayer outer membrane forms, of which the individual layers are exfoliated by nucleating CC crystals. 136,137 Furthermore, recent experimental works show a thin ACC layer on the surface of the CCC tablets of the Figure 11. Growth rates predicted by models PF1 (red lines) and PF2 (blue lines) for the ion-by-ion mechanism from aqueous solution and for crystallization from ACC precursor. Note the ∼7 orders of magnitude difference in growth velocity predicted for the two mechanisms. For comparison, the range of experimental data is also shown (gray domain). Note the oscillating growth rate in the layered domain.

JACS Au
pubs.acs.org/jacsau Article nacre, 141,142 implying that (h)ACC particles may play an essential role in the formation of nacre. This, in turn, indicates that a particle-by-particle addition is more appropriate for solidification, as in the case of coral skeletons (see discussion in section 3.3). Summarizing, while models PF1 and PF2 cannot capture all details of the formation mechanism of nacre, remarkably similar microstructures are generated. This raises the possibility that on the basis of the mechanisms these models realize (diffusion controlled solidification at high driving forces), one may design and prepare artificial mollusk shell structures that inherit the mechanical excellence emerging from the hierarchical sequence of the granular, prismatic, and nacre-like ultrastructures. Thus, the present work opens up the way toward a novel design strategy for creating biomimetic/bioinspired composite materials.
Finally, to address the evolution of columnar nacre within the phase-field theory, one can incorporate preexisting organic membranes into the simulations "by hand", including the mineral bridges/aligned holes seen in the experiments. 33,34,138,140 We used PF1 to model the formation of CCC stacks shown in Figure 12. The thin organic walls were assumed to have amorphous structure (local orientation varied randomly pixelwise) and we applied a boundary condition that ensured a contact angle of 100°with the solid−liquid interface. Simulations of this kind yield 2D "pyramid-like" stacks of CC layers as shown in Figure 13. Such simulations can be used to explore the effect of such parameters as growth anisotropy, contact angle, hole size/position, and so on. For example, Figure  13b implies that the growth velocity of stack height increases with increasing hole width.

Helical Structures Predicted by Model PF3
Spectacular screw dislocation-like helical structures have been observed in mollusk shells akin to patterns formed in oscillating chemical reactions (Figure 4). These 3D structures cannot be addressed in models PF1 and PF2 as the scalar orientation field used in them is valid in only 2D (as in 3D minimum the three Euler angle are needed to define the crystallographic orientation). Therefore, we use model PF3 to explore the possibility of forming such objects within the framework of the phase-field theory. Here two solid phases (α and β) precipitate from a homogeneous ternary liquid. Owing to the lack of relevant information, we retain the materials parameters used in ref 102. Under appropriate conditions, shown by the green diamonds in Figure 14, a layerwise structure composed of alternating α (mineral) and β (organic) layers form (the third component is water, the crystal grows into hACC), an analogue of the "nacre" observed in model PF2. Note that the layer-bylayer growth mechanism, by which the layered structure forms is spinodal nucleation of one phase on the other, and thus can be viewed as an extreme case of island growth. In this regime, we see the formation of helical structures that emerge in pairs of opposite chirality. As the spiraling eutectic dendrites in ref 102, these defects originate from an instability associated with diffusion of the third component.
Different views and sections of such structures are shown in Figure 15. This chiral structure closely resembles the screw dislocation-like defects reported in experiments. 14,39 The similarity could be enhanced by incorporating kinetic/interfacial Here ss and s stand for "surface sheets" and "organic sheet", whereas sm and ilm indicate "surface membrane" and "interlamellar membrane", respectively.

JACS Au
pubs.acs.org/jacsau Article energy anisotropies yielding faceted growth perimeters. It cannot, however, be excluded that the screw dislocation-like helical structures have here a different origin than in the reality. For example, our dislocation pairs do not recombine even in the long time limit, presumably because of the lack of mechanical stresses that are not incorporated into the phase-field models used in this study. Further investigations are yet needed to clarify whether phase-field models of island growth 143 or anisotropic pinning 144 could provide more realistic dynamics. Simultaneous formation of alternating organic and mineral layers in the vicinity of topological defects predicted by our simulations is a working hypothesis. Previous work suggests that the formation of spiraling organic membranes precedes mineral deposition. 39,145 Biological systems often use liquid crystals based on chitin as a template. 146 It is hypothesized 39,40,145 in this case that self-organization of the liquid crystal phase leads to the formation of a helical organic scaffold, which serves as a template for the mineralized structures of the same configuration. This mechanism of scaffold formation may be analogous to the formation of helical structures in Liesegang systems, 147,148 which raises the possibility of using the phase-field inventory developed for such systems. 149 Note that the present simulations were not optimized for the case of nacre. Further work is underway to characterize these structures and the dynamics of their formation within the phasefield theory.

Modeling of Coral Skeletons in Model PF1
Next, using model PF1, we try to find a qualitative answer to the question of why the skeleton of some corals species contain small randomly oriented crystallites, "sprinkles" (see Figure 4e), which occur at the perimeter and along grain boundaries and even form bands, whereas other species do not display this behavior, an observation discussed in some detail in ref 28. In our previous work, we demonstrated that conditions of mineralization can influence the amount of sprinkles; however, we have addressed only tangentially the formation of sprinkle bands.
In this section, we address the latter phenomenon. We hypothesize that the coral polyp sits on the corallite (top of the skeleton) and emits a supersaturated extracellular calcifying fluid (CF), which is not in direct contact with the seawater but fills in the cavities of the porous skeleton. Recent work indicates that the coral skeletons are deposited biologically and actively via attachment of hACC nanoparticles (of diameter 50 to 400 nm), while ion-by-ion addition fills the interstitial space among them. 150 This leads to the formation of a thin (<1 μm)  As the hACC/ACC surface layer remains thin, crystallization is apparently fast enough to keep pace with hACC deposition. Thus, diffusion of the hACC nanoparticles is expected to be the rate-limiting phenomenon. Thus, although particle diffusion does not influence the thermodynamic driving force of the ACC → CCC transition, it does influence the effective crystal growth rate. The effective rate is small on the sides of the fingers due to the low concentration of nanoparticles and allows a longer time the consolidation of the crystal layer yielding larger grain sizes. At the tips, the effective crystal growth rate is larger because of a higher concentration of the nanoparticles; thus, a shorter time is available for consolidation, and the crystal grains remain smaller. This is observed in the experimental images (Figure 16a−c).
The diffusion coefficient for hACC particles of radius R p = 150 nm can be estimated by the Stokes−Einstein relationship D p = k B T/(6πR p η), which leads to D p ≈ 2 × 10 −12 m 2 /s, for T = 300 K and a viscosity of η = 1 mPa·s. Taking a medium growth rate 61−63 of v = 5 cm/yr ≈ 1.6 × 10 −10 m/s for the coral skeleton, one obtains a diffusion length of l D ≈ 0.92 mm, whereas for fastgrowing corals (37 cm/yr) 63 l D ≈ 0.12 mm, which do not rule out diffusion-controlled solidification.
Accordingly, model PF1 is applied here as follows: the phasefield and the coupled (particle-) concentration field control solidification on the time scale of particle diffusion and lead to a Mullins−Sekerka-type diffusional instability 151 of the interface. The magnitude of the orientational mobility, in turn, determines whether the forming solid becomes orientationally homogeneous (single crystal), partly ordered (polycrystalline), or fully disordered (amorphous).
Since in the case of coral skeletons we do not need to produce nacre-like alternating organic−inorganic layers via a transition from diffusionless to diffusive growth, in model PF1 we can use conditions, where dv/D c,M ≪ 1; that is, the system is far from the diffusionless growth regime. The dimensionless model parameters and boundary conditions are the same as those we used in ref 28 for coral skeletons.
The polycrystalline structure emerging under these conditions due to the diffusional instability is shown in Figure 16d− f. Note the rough surface and the liquid channels, with small crystallites at the tips and along the spine of the branches, and larger crystallites at the lateral surfaces. This distribution of grain size is the result of the fact that the growth velocity at the tips is larger because of the larger supersaturation (the tip meets fresh, nondepleted liquid, see Figure 16f), whereas in the interarm channels, the fluid is depleted, so growth is slow, yielding larger crystallites. This phenomenon is the result of growth front nucleation (GFN) that leads to more frequent GFN events with increasing growth velocity as discussed in refs 75−78,152. Whether the combination of diffusional instability with GFN or some other biology directed mechanisms is responsible for the appearance of the rough surface of the skeleton is unclear at present. It is thus desirable to investigate further consequences of the hypothesized control of the grain size distribution via diffusional instabilities and GFN.
Along these lines, we make predictions using model PF1 that may be tested experimentally and validate or disprove our assumptions. For this reason, we investigate the effect of model parameters that influence the amount of sprinkles forming in the model and thus may offer explanation why it varies from one species to the other.
Previous work indicated that the orientational mobility (related to the rotational diffusion coefficient), the thermodynamic driving force (influenced by supersaturation or temperature) may influence the intensity of GFN (formation of new grains at the growth front). [75][76][77][78]153 In the framework of this study, we varied individually these parameters, and in all cases, we obtained a transition from microstructures dominated by sprinkles to microstructures with few or no sprinkles. Of these parameters, only the temperature can be controlled with relative ease. Variation of particle concentration in the calcifying fluid below the coral polyp, or controlling the rotational diffusion coefficient of the particles in the fluid are probably beyond the reach of the experimenter. Therefore, we present only the microstructural/morphological changes predicted as a function of temperature (see Figure 17). Apparently, according to our model, the coral skeletons grown at low temperatures should have larger amount of sprinkles than those grown at higher temperatures. This finding raises the question whether the amount of sprinkles is indeed a characteristic feature of the individual coral species (i.e., determined biologically) or some differences in the circumstances of mineralization (temperature and/or supersaturation) are responsible for the deviations. In

JACS Au
pubs.acs.org/jacsau Article any event, these simulations may offer a natural explanation for the differences seen in the amount of sprinkles in the skeleton of different coral species.

CONCLUSIONS
We have demonstrated that qualitative coarse-grained phasefield modeling offers a methodology to address specific mesoscale aspects of nonclassical crystallization phenomena taking place during biomineralization. In particular, we investigated to what extent qualitative phase-field modeling can contribute to the understanding of microstructure evolution during the formation of mollusk shells and coral skeletons for various species. We applied three different phase-field models: PF1, PF2, and PF3. Our present findings can be summarized as follows: (1) Ultrastructure specif ic to the shells of mollusks Unio pictorum, Nautilus pompilius, and Haliotis asinina: Driving the solidification process from solute trapping toward partitioning via decreasing the thermodynamic driving force, binary phase-field models PF1 and PF2 recover the common sequence of granular → prismatic → nacre ultrastructures on a reasonable time scale, if CC crystallization takes place via an amorphous precursor.
In contrast, within this scenario, CC crystallization via ion-by-ion deposition from aqueous solution appears to be orders of magnitude too fast when compared to experiments. (2) Nacre formation in mollusk shells: Models PF1 and PF2 describe reasonably well the formation of not only the granular and prismatic domains, but the appearance of sheet nacre as well, in which case the models indicate alternating precipitation of the organic and mineral components. The models seem to reproduce even such details as mineral bridges and aligned holes. Yet, for obvious reasons, they cannot predict the formation mechanism of columnar nacre, in which the formation of organic membranes precedes CC precipitation. However, representing the preexisting organic membranes via appropriate boundary conditions, a reasonable description can be obtained even in this case. (3) Screw dislocations in mollusk shells: Ternary phase field model PF3 predicts the formation of screw dislocations pairs in 3D, a phenomenon analogous to the experimental findings. Inclusion of elasticity into the model is needed to capture the proper dynamic behavior during growth. (4) Sprinkle formation in coral skeletons: model PF1 was used to explore the possible mechanism for the formation of nanoscale crystallites "sprinkles", whose presence was reported recently in the skeletons of certain coral species. Assuming a diffusion controlled mechanism in confined space, we observe the formation of sprinkle bands at the spine of the arms of the corallite as a trace of fast solidification at the arm tips, whereas larger crystallites form at the sides of the arms. The simulations show that varying the orientation mobility (proportional to the rotational diffusion coefficient of the molecules/ions) or the driving force of crystallization (via changing either the supersaturation or the temperature), one can control the amount of sprinkles between essentially no sprinkle and dominantly sprinkled microstructures. Unquestionably, the applied models should be viewed as only minimum models of the processes taking place during biomineralization. Yet, we believe that phase-field modeling complemented with biochemical/biological information has the potential to contribute to a better qualitative or even quantitative understanding of morphogenesis in simple cases. Introduction of more complex models can certainly improve the mathematical representation of the associated phenomena. Ultimate limitations of such approaches stem from the fact that living organisms cannot be modeled within this framework: they can only be represented by boundary conditions of different complexity. Despite these, in specific cases, we recovered structures closely resembling their biogenic counterpart. The resemblance of the simulated and natural biominerals suggests that, underneath the immense biological complexity observed in living organisms, the underlying design principles for biological structures may be so simple that they can be understood with simple math and simulated by phase-field theory.
Finally, we note that our simulations outline conditions, under which standard materials science processes can be used to create inorganic substances that mimic the microstructures observed to form in living organisms, a knowledge that may open up ways for creating new biomimetic/bioinspired composite materials.