Spontaneous mirror symmetry breaking in heterocatalytically coupled enantioselective replicators† †Electronic supplementary information (ESI) available: Methods and text on stability analysis. Five figures showing simulation examples. See DOI: 10.1039/c6sc02446g Click here for additional data file.

Hypercycles proposed as a chemical basis for the selection of biological replicators may lead to homochirality when fed from achiral resources.


Methods
Simulations were performed by numerical integration of the differential rate equations according to rate-equation theory as applied to chemical kinetics (mean field assumption). The concentration units are mol L -1 and the different rate constants have the rate values in units of mol s -1 . Numerical integration was performed with the Mathematica® program package. For a set of parameters corresponding to the system at, or very near to, the bifurcation point the numerical integration is highly sensitive to minute differences between the reaction parameters, so that the inherent numerical noise of the calculations suffices to bifurcate the system towards a chiral outcome or made insensitive to SMSB. In our simulations we have suppressed this computational noise, arising from round-off errors, by setting a high numerical precision of the input parameters (100 significant decimal digits and exact number representation of the reaction rates and the initial concentration values (for example "1 + 1·x 10 -2 " instead of "1.01" or "1. + 1 x 10 -10 " or "1 + 1. x 10 -2 "). Integration methods of "StiffnessSwitching" and "WorkingPrecision" of up to 50 were used in the present calculations.
The ee is expressed in % by 100 x (([L] - [D])/ (([L] + [D])). The fluctuations of chirality able to take the system out from the racemic branch were simulated by using an initial ee of products/catalysts lower than that expected from the statistical fluctuations about the ideal racemic composition, i.e an initial ee (%) < 67.43 x (N -0.5 ), where N is the number of chiral molecule. S1 In the examples shown here the initial ee was an order of magnitude lower than that expected from the chiral statistical fluctuation. The numerical integration was run between 0 s to 1 x 10 20 s. This limit of time, three orders of magnitude larger than the age of the universe, allows us to estimate the concentration value of the stationary state. In the simulations presented here SMSB occur between one day (of chemical significance in potential applied synthesis) and a few decades of years. The effect of an external chiral permanent perturbation of the transition states can be simulated taking k iD ≠ k iL . In this case values as low as k iD /k iL = 1 x 10 -9 (transition free energy differences at room temperature in the order of mJ mol L -1 ) suffice to drive the systems at SMSB conditions from the racemic to the chiral state.
The search of the reaction parameters leading to SMSB was, in the case of the two-replicator system of Figs. 2, S3 and S4, previously estimated from the algebraic adimensional analysis conclusions (Fig. S2.). However, for the compartmentalized system at two temperatures and for hypercycle composed by several replicators the parameter region for SMSB was searched by trial and error of numerical simulations. A complex behavior was observed, for example in the chemical mass, that in simple Frank like systems represents a value below or above which SMSB can be obtained. In the systems studied here, we found regions of chiral solutions bounded by two critical values of the total chemical mass.

Stability studies 1 Single Chiral Replicator: Quadratic autocatalysis
We consider direct production of the replicator enantiomers (R L and R D ) from an achiral source : [S1] k 0 and quadratic autocatalysis driven by constant concentration external reagents X,Y: such that the number of "internal" molecules A+R L +R D is conserved. These correspond to the transformations [1] and [2] of Scheme 1 in the main text. The temporal and T spatial physical dimensions of the reaction rate constants (indicated here by the V brackets) [S3] and [S4] allow us to express the associated kinetic rate equations in terms of dimensionless rates and concentrations [1,2]. Changing variables to  = R D +R L , [S1,S2] then imply [see for example refs. S2 and S3 for more details]. [S5] [S6] where  = k 0 t is dimensionless time, the dimensionless concentrations, and [S7] Equations (S5,S6) admit one racemic and two chiral stationary solutions: [S8] [S9] To assess dynamic stability, we linearize the rate equations [S5,S6] and consider arbitrary fluctuations around the stationary solutions. Their time dependence is determined by [S10] where [S11] We evaluate the array over one of the stationary solutions [S8,S9] (racemic or chiral) A and compute the corresponding pair ( 1 ,  2 ) of eigenvalues. Doing so for the racemic solution we find: [S12] whereas for both the chiral solutions we find [S13] This establishes that quadratic autocatalysis together with direct production yields the racemic solution as the unique final stable state, despite the fact that the autocatalysis is driven by external reagents.

Single Chiral Replicator: Cubic autocatalysis
We consider direct production of the replicator enantiomers (R L and R D ) from an achiral source : [S14] k0 and quadratic autocatalysis driven by constant concentration external reagents X,Y: [S15] k R such that the number of "internal" molecules A+R L +R D is conserved. These correspond to the transformations [1] and [2] of Scheme 1 in the main text. The temporal and T spatial physical dimensions of the reaction rate constants (indicated here by the V brackets) [S16 and [S17] allow us to express the associated kinetic rate equations in terms of dimensionless rates and concentrations [1,2]. Changing variables to  = R D +R L , [S14,S15] then imply [see for example refs. 25,26 for more details]. [S18] [S19] where  = k 0 t is dimensionless time and [S20] The dimensionless concentrations appearing in [S18,S19] are obtained by rescaling their dimensionfull counterparts as follows .

No racemization no SMSB
 If we omit the forward/reverse direct production steps [S14], then proceed with the non-dimensionalization, the stationary solution of the resultant [S18] is for y = 0 [and for which [S19] holds identically). The corresponding eigenvalues for time-depende fluctuations  and y around this racemic state are Both   < 0,  y are negative when evaluated on the above solution for , implying the  stability of the mirror symmetric state. There is thus no SMSB in the absence of the direct production steps, regardless of the values of and the total system concentration . g C

Racemization: Phase space for SMSB
When the direct production steps are reinstated, when y = 0 we must solve the cubic equation [S23] for the positive solutions  > 0.
The pair of eigenvalues are given by [S24] [S25] Now can be positive on the positive solutions to [S23]. The phase space for SMSB is y  described by the three parameters u, g, C defined in [S20]. A representative region where mirror symmetry is broken ( y > 0) is depicted in the phase diagram Fig. S1. Figure S1. Effect of varying total system concentration C, the ratio of the rates of u inverse to direct production and the ratio of the rates of inverse to forward cubic g autocatalysis (driven by external reagents). Left: the three dimensional region represents points (colored) for which the racemic state y = 0 is unstable to small perturbations (y >0) and hence where SMSB results. Right: a transverse "slice" of this region at C =2.5 indicting the nonlinear dependence of on and the maximum value of the former. g u Left: That maximum increases with . For a given system concentration C > C min , there C is a finite range in both and for which SMSB occurs. Now u  0 corresponds a u g rapid racemization or direct production rate with respect to its inverse: k 0 >> k -0 . The finite upper bound on corresponds to the opposite situation, namely when k -0 << k 0 But u the forward rate of direct production cannot be too small, otherwise there is no SMSB. 0 k The limit g  0 corresponds to relative fast rates of the forward to reverse (driven) cubic catalysis. For a given , there is a finite upper bound on indicating minimum values of u g (driven) forward catalysis needed to set off SMSB. Note this range of minimal rates of forward (driven) catalysis is maximized for intermediate values of (right hand side u Fig.S1.). Increasing widens both the range in and the maximum value of for C u g which SMSB can occur. where dimensionless time-parameters are  = k 0j t for j = 1,2, and

The Two Replicator Hypercycle: replicators aided by homochiral cross-catalysis
All the concentrations in the above rate equations are dimensionless. The relation between the dimensionless and dimensionfull concentrations (e.g., [A]) is given by [S35] The total system concentration is conserved: C [S36] The dynamics is thus described by the six independent dimensionless parameters in [S34,S36].
We seek stationary solutions (1),  (2) corresponding to the global racemic state y (1) = y (2) = 0. In this configuration, the fluctuations in the net chiral mass in each species  (1) ,  (2) decouple from the chiral fluctuations about the racemic and so the (2) (1) , y y   4x4 Jacobian matrix decomposes into block-diagonal form with the two 2x2 sub-blocks: The time dependence of the matter and chiral fluctuations about the global racemic obeys The signs of the eigenvalues of the 2x2 arrays and evaluated over these acceptable A B solutions indicate whether the racemic solution is stable or not. If any one of the (four) eigenvalues is positive, the racemic state is unstable to arbitrary small perturbations which thus trigger the onset of SMSB. A parameter survey leads to the representative phase diagram discussed below.

Phase Space for SMSB:  = u, h = g, s = 1 and C
The full phase space of this two-replicator system is six-dimensional, nevertheless, we can obtain considerable insight into the roles played by the total system concentration, the racemization and the rates of the (driven) heterocatalyses by specializing to a three dimensional parameter space. To this end, we here set v = u, h = g, s =1 and then vary u, g, C. Note that corresponds to the ratio of the rates of inverse to h forward catalysis for the species R (2) when s =1, and in this limit is analogous to the parameter for species R (2) : g An illustrative region where mirror symmetry is broken is depicted in the phase diagram Fig. S2. This is qualitatively similar to the phase portrait for the single replicator case discussed above Fig. S1. A notable quantitative distinction is that the minimum value of required for SMSB or two-replicators has increased with respect to the single C replicator case. Figure S2. Effect of varying total system concentration C, versus the ratio u =  of the rates of inverse to direct production of both species R (1) , R (2) and the ratio g = h of the rates of reverse to forward autocatalysis of both species (driven by the external reagents). Here s =1. The three dimensional figure represents points for which the global racemic state y (1) = y (2) = 0 is unstable to small time-dependent perturbations. Figure S3. SMSB in the system of Fig. 1 but in a reaction network including also the non-cross catalyzed autocatalysis (2) and the same achiral resource A for the two chiral replicators ( 1 R D ( 1 R L ) and 2 R D ( 2 R L )). The autocatalysis (2) as well as (3) are driven by the reagent X/Y at constant concentrations. (1 x 10 -3 + 1 x 10 -23 ) mol L -1 (see Methods) The initial concentration of the rest of the species was set to zero. The same final state when starting from a different relationship of initial concentrations, i.e. the final species composition corresponds to a thermodynamically controlled stationary state. Figure S5. SMSB of a two replicator chiral hypercycle in a closed system with permanent nonhomogeneous temperature distribution: two compartments at different temperature exchanging the same volume of solution. For a more detailed decription of this type of systems see for details Astrobiology 2013, 13, 132-142. The change on equilibrium constants would agree with that expected for temperature differences as those of deep ocean hydrothermal wells (411° K and 275° K) and for an exothermic transformation showing negative entropy, as expected for replicator obtained by polymerizating reactions.  The ee output value (given as L enantiomers excess value) is that of full homochirality. Reaction rate constants (k iR /k -iR ): 1 R: 1 x 10 4 /10; 2R: 1 x 10 3 /0.5. Initial resource concentrations in the reactor and in the constant input volume (0.1 µL s-1): [A] o = 1 x 10 -4 mol L -1 ; [B] o = 1.5 x 10 -4 mol L -1 . Initial replicator concentrations in the reactor were 1 x 10 -6 mol L-1 and the initial chiral fluctuation was simulated by an additional concentration of 1x10 -23 mol L -1 in the L enantiomers.