Equilibrium phases and domain growth kinetics of calamitic liquid crystals

Prior to investigating the phase ordering kinetics, we first visit the problem of equilibrium phase transitions in the GB model to precisely identify the quench temperatures for the Nm and the Sm regimes. To sample the available phase space, MC simulations are performed in the canonical ($NVT$) ensemble. We use the simulation program DL_MONTE Purton et al. (2013); Brukhno et al. (2021) for this purpose. DL_MONTE is a general-purpose MC program which supports a wide range of MC simulation techniques and interatomic potentials. However, prior to this work it was not applicable to particles with implicit orientations and anisotropic interaction potentials, for instance the GB particles. Our work is a step in that direction since it involved extending DL_MONTE to make it suitable for such systems. The latest version of DL_MONTE which includes our improvements is available at DL_ . This code is a beneficial resource for the research community to study uniaxial GB systems with a number of MC simulation techniques, including grand-canonical and constant-$NPT$ ensembles Frenkel and Smit (2002). In this work, we only present results using MC in the $NVT$ ensemble. Our DL_MONTE input files with all the necessary commands, parameters and comments are provided in the Supplementary Material for interested readers together with further details regarding the improvements we have made to DL_MONTE to treat the GB and similar models.

We consider a system of $N$ ellipsoidal particles interacting via the GB potential specified in Eq. (1) and use $\kappa=3.0$, $\kappa^{\prime}=5.0$, $\mu=1.0$ and $\nu=3.0$ unless specified. It is convenient to define scaled variables $E^{\ast}=E/\epsilon_{0}=\sum_{i,j;i<j}^{N}(E_{ij}/\epsilon_{0}$), $E_{1}^{\ast}=E^{\ast}/N$, $T^{\ast}=k_{B}T/\epsilon_{0}$, $r^{\ast}=r/\sigma_{0}$ and $\rho^{\ast}=N\sigma_{0}^{3}/V$. Simulations have been performed on cubic lattices with $N=512$ and $1000$ particles using periodic boundary conditions and $\rho^{\ast}$ is set to 0.30. The initial configuration is a perfectly aligned state. The MC moves are performed according to the Metropolis algorithm Frenkel and Smit (2002). In these moves, a randomly selected LC molecule is either translated or rotated with translations and rotations attempted in each move with equal probability. The maximum angle to rotate an LC molecule is kept as $15^{\circ}$. In order to compute the energy efficiently, a spherical cutoff of $r_{c}=4.0\sigma_{0}$ is employed in conjunction with the Verlet neighbor list scheme Frenkel and Smit (2002). The simulations are performed for $10^{6}$ MC cycles (one MC cycle corresponds to $N$ MC moves). The initial $5\times 10^{5}$ cycles are necessary for equilibration. The remaining $5\times 10^{5}$ cycles are used for thermal (block) averaging of various thermodynamic quantities of interest and the ensemble averaging is performed over 500 configurations.

The GB phase diagram is obtained by studying the average energy per particle $\langle E_{1}^{\ast}\rangle$ vs. $T^{\ast}$ for values of the energy anisotropy $\kappa^{\prime}$ = 1.25, 2.5, 5.0, 10.0 and 20.0. For each value, the transition temperatures $T_{c}^{1}$ (I$\rightarrow$Nm transition) and $T_{c}^{2}$ (Nm$\rightarrow$Sm transition) are identified from the discontinuity in $\langle E_{1}^{\ast}\rangle$. It is pertinent to point out here that for our coarsening experiments, we only need approximate boundaries as we consider the quenches far away from these. The I and Nm phases are confirmed by evaluation of the orientational order parameter, $\mathcal{S}$ given by:

where $\cos\theta_{i}=\bm{u}_{i}\cdot\bm{n}$ and the angular brackets $\langle\cdot\cdot\cdot\rangle$ indicate an ensemble average. $\mathcal{S}=0$ in the I phase, while $\mathcal{S}=1$ in the perfectly aligned Nm phase. The defects correspond to regions with $\mathcal{S}\simeq 0$, even if the defect cores are biaxial Callan-Jones et al. (2006).

The smectic (layered) phases can be distinguished from the I and Nm phases by evaluating the translational order parameter Bates and Luckhurst (1999):

where $\bm{r_{i}}$ is the position vector of the LC molecule $i$ and $d_{l}$ is the layer spacing. In simulations, $\mathcal{T}$ is determined by first separately performing the ensemble averaging of the real and imaginary terms: $\cos(2i\pi\bm{r_{i}}\cdot\bm{n}/d_{l})$ and $\sin(2i\pi\bm{r_{i}}\cdot\bm{n}/d_{l})$, followed by calculation of the modulus for the ensemble average and further maximizing it with respect to $d_{l}$ Bates and Luckhurst (1999). In the I phase, $\mathcal{T}$ $\rightarrow$ 0 as the layers are not well defined. On the other hand, $\mathcal{T}$ $\rightarrow$ 1 in a perfectly layered structure.

A clear distinction between the SmA and SmB phases can be made by evaluating the hexatic bond-orientational order parameter given by Bates and Luckhurst (1999):

where the summation is over the nearest neighbours (nn) $k$ of the LC molecule $i$, $\phi_{ik}$ is the angle between the vector ($\bm{r_{i}}-\bm{r_{k}}$) projected onto the plane normal to the director and a fixed reference axis ($x$ axis say) and $w(r_{ik}^{\ast})$ is a cutoff function to select the nn for evaluation of $\psi_{6}(\bm{r_{i}})$. (Like $\mathcal{T}$, for $\mathcal{C}_{6}$ also first the ensemble averaging is done separately for the real and imaginary terms followed by calculation of the modulus for the average.) It is important to use the cutoff function as the number of nn might not be 6 and could be 7, 5 or 4 when the local translational order is imperfect. In our work, we used the procedure in reference Bates and Luckhurst (1999) to evaluate this cutoff function: $w(r_{ik}^{\ast})$ is unity for $r_{ik}^{\ast}$ below 1.4, zero for $r_{ik}^{\ast}$ above 1.8 and with a linear interpolation in between these two extremes. If there is hexatic order, $\mathcal{C}_{6}$ has an appreciable non-zero value and also the ensemble average for the cutoff function, $\left\langle w(r_{ik}^{\ast})\right\rangle\rightarrow 6$. In it’s absence, $\mathcal{C}_{6}$ vanishes and $\left\langle w(r_{ik}^{\ast})\right\rangle<6$.

The radial distribution function $g(r)$ is routinely obtained in scattering experiments, and measures the probability of finding two molecules separated by distance $r$ relative to that in an ideal gas Frenkel and Smit (2002). It is a useful tool to distinguish between the local order in the different phases Bates and Luckhurst (1999), and is given by:

where $\rho_{0}=N/V$ is the density of the ideal gas and $\overline{\rho(r^{\ast})}$ is the average density of the system around $r^{\ast}$. The numerical evaluation is facilitated by the following formula Frenkel and Smit (2002); Billeter and Pelcovits (1998):

The $\delta$ function is unity if $r_{ij}^{\ast}$ falls within the shell centered on $r^{\ast}$ and is zero otherwise. The division by $N$ is done to normalize $g(r^{\ast})$ to a per-molecule function. By construction, $g(r^{\ast})=1$ for an ideal gas and any deviation implies correlations between the particles due to the intermolecular interactions. In the Nm phase, it has a noticeable maximum for small-$r^{\ast}$ and is 1 for large-$r^{\ast}$ indicating short-range order. In the SmA phase, $g(r^{\ast})$ exhibits a large nn peak and small oscillations thereafter with periodicity $\sim d_{l}$ due to the layering. In the SmB-H phase on the other hand, $g(r^{\ast})$ shows an oscillatory behaviour with many sharp peaks separated by $r^{*}\simeq 1$ and a split-peak at $r^{*}\simeq 2$ characteristic of hexatic order within the layers. In the case of SmB-C order, the peak intensities do not decay till large distances.

We now proceed to evaluate the phase diagram in $(T^{\ast},\kappa^{\prime})$-space for $\kappa^{\prime}$ = 1.25, 2.5, 5.0, 10.0 and 20.0. In Fig. 2(a), we show the variation of the per-particle equilibrium energy $\langle E_{1}^{\ast}\rangle$ vs. $T^{\ast}$ for $N=512$ (open up triangles) and $N=1000$ (open down triangles) at $\kappa^{\prime}=5.0$. (For brevity, we do not show the data sets for other values of $\kappa^{\prime}$.) The angular brackets indicate thermal averages. The simulation data coincide, indicating that the average energies are independent of the system size. We also plot the corresponding benchmarking results of Berardi et al. for $N=512$ (white circles) and $N=1000$ during heating (white up triangles) and cooling (white down triangles) protocols Berardi et al. (1993). There is an excellent agreement with our simulation results obtained using DL_MONTE, even though the initial conditions in the two sets of simulation experiments are quite distinct. Furthermore, the effects of system size on the average energy are negligible (except near the transition in some cases). The left and right edges of the Nm region (in green) provide the scaled transition temperatures $T^{\ast 2}_{c}$ and $T^{\ast 1}_{c}$ as indicated. Both the transitions are discontinuous (first-order).

To confirm the phase corresponding to each data point, we have also evaluated for the $N=1000$ system, the temperature variation of the average order parameters $\langle\mathcal{S}\rangle$, $\langle\mathcal{T}\rangle$ and $\langle\mathcal{C}_{6}\rangle$, or $g(r^{\ast})$ vs. $r^{\ast}$ as appropriate. Fig. 2(b) shows the typical behaviour of the order parameters in various phases. Note the rise in $\langle\mathcal{S}\rangle$ as the temperature crosses $T^{\ast 1}_{c}$ from above, indicating the onset of nematic order. Similarly, note the rise in translational order parameter $\langle\mathcal{T}\rangle$ as the temperature crosses $T^{\ast 2}_{c}$ indicating the development of SmA order (while evaluating $\mathcal{T}$, we have observed that the layer spacing $d_{l}$ which maximizes $\mathcal{T}$ is less than the aspect-ratio $\kappa=3.0$ due to interdigitation of the layers Bates and Luckhurst (1999)). At slightly lower values of $T$, the hexatic order parameter $\langle\mathcal{C}_{6}\rangle$ becomes significant suggesting that the phase is SmB-H. Accurate determination of the phase boundaries will require careful evaluations. We refrain from going in this direction as the focus of the present study is on coarsening and approximate phase boundaries are sufficient for this purpose. (We wish to mention here that though the equilibrated configurations, their energies in various phases and the order parameter $\langle\mathcal{S}\rangle$ were obtained using DL_MONTE, the codes for evaluation of the order parameters $\mathcal{T}$ and $\mathcal{C}_{6}$ were written separately.)

To evaluate $g(r^{\ast})$, we have used a set of $745$ bins with a separation cut-off of $r^{\ast}=7.45$ which yields $\Delta r^{\ast}=0.01$ (see Supplementary Material for details regarding specification of parameters for DL_MONTE). Fig. 2(c) shows the $g(r^{\ast})$ vs. $r^{\ast}$ behaviour for the Nm phase (green line, $T^{\ast}=2.5$) and the Sm phase (blue line, $T^{\ast}=1.0$). It is characterized by multiple sharp peaks with a split peak at $r^{\ast}\simeq 2$, characteristic of the SmB-H order. Further, it can be seen that $g(r^{\ast})$ decays as $r^{\ast}$ increases, implying that the translational order is lost at larger distances and hence the SmB phase is not of crystal type.

Fig. 2(d) depicts the variation of $T^{\ast 1}_{c}$ (white up triangles) and $T^{\ast 2}_{c}$ (white down triangles) for a range of $\kappa^{\prime}$ values. To be noted here is that $T^{\ast 1}_{c}$ decreases considerably with increasing $\kappa^{\prime}$ thereby shrinking the Nm phase. On the other hand, $T^{\ast 2}_{c}$ increases only slightly. With the original GB parameterization, de Miguel et al. also reported similar observations de Miguel et al. (1996).

We now turn to a study of domain growth kinetics in the $d=3$ GB model via MD simulations. To the best of our knowledge, they have not been addressed. There have been some coarsening studies for the Nm phase of the GB model in $d=3$. A relevant contribution in this context is by Billeter et al. Billeter et al. (1999). They observed the usual defect structures (e.g. disclinations and monopoles), but could not extract a reliable growth law due to the small system size considered ($\sim$60000 particles). Using large systems ($\sim$260000 particles), we perform deep quenches from the I phase to the Nm and the Sm phases [shown in Fig. 2(d)] and allow the system to evolve for long times. Our main interest is to determine the novel defects and growth laws in the Sm phase, and the impact of the energy anisotropy $\kappa^{\prime}$ on ordering. We also undertake analogous studies for the Nm phase to complement previous works Billeter et al. (1999).

All our MD simulations have been performed in the $NVT$ ensemble using the LAMMPS software package Plimpton (1995); Lam . The details regarding implementation of the model into LAMMPS ¹¹1Note that LAMMPS technically implements the generalized Gay-Berne model Berardi_1995 in which interactions of the ellipsoidal particles are characterized via a (diagonal) shape matrix $S=$ diag($\sigma_{a},\sigma_{b},\sigma_{c}$) and a (diagonal) energy matrix $E=$ diag($\epsilon_{a},\epsilon_{b},\epsilon_{c}$), where $\sigma_{a},\sigma_{b},\sigma_{c}$ are the lengths and $\epsilon_{a},\epsilon_{b},\epsilon_{c}$ are the relative well depths of interaction along the three semi-axes of an ellipsoid. The GB model discussed in Sec. II and employed in Sec. III was the potential originally presented for ellipsoidal particles of equal size. The generalized GB model was developed to represent dissimilar biaxial ellipsoids Berardi_1995. It reduces to the GB model described in Sec. II when the molecules become uniaxial. Hence, for our study, we work with the generalized GB model in the uniaxial limit. and also the analytical expressions for the forces and torques, have been described in Brown et al. (2009). We consider $N=262144$ uniaxial ellipsoidal particles confined in a cubic box of linear size $L_{s}\sigma_{0}$ with periodic boundary conditions in all three coordinate directions. Hence the volume $V=$ ($L_{s}\sigma_{0}$)³ = ($95.6033\sigma_{0}$)³ such that $\rho^{\ast}\simeq 0.3$. The parameters in the GB potential of Eq. (1) are those employed by Berardi et al. Berardi et al. (1993). We also study the effects of varying the energy anisotropy parameter $\kappa^{\prime}$. The MD runs are carried out using the standard velocity Verlet algorithm. In LAMMPS, the dimensionless MD time unit $t_{0}=\sqrt{m_{0}\sigma_{0}^{2}/\epsilon_{0}}=1.0$. We choose the reduced MD integration time step $\Delta t^{\ast}=\Delta t/t_{0}=0.001$. The temperature $T^{\ast}$ is controlled and maintained constant via the Nosé-Hoover thermostat, which is known to preserve hydrodynamics Frenkel and Smit (2002); Nosé (1984); Binder et al. (2004). The homogeneous initial configurations are prepared by equilibrating the system at a high temperature ($T=6.0$) for about $10^{5}$ MD steps. To initiate the coarsening process ($t=0$), the system is quenched to the indicated temperatures in Fig. 2(b). The evolution of the system is then monitored. All statistical quantities of interest are averaged over $19$ independent initial conditions. Our input files and parameters used in LAMMPS are provided in the Supplementary Material.

For a translationally invariant system, the usual probe to characterize configurational morphologies is the equal-time correlation function Puri and Wadhawan (2009):

where $\psi(\vec{r},t)$ is a suitable order parameter, $\vec{r}=\vec{r_{2}}-\vec{r_{1}}$ and $\langle\cdot\cdot\cdot\rangle$ represents the ensemble average. Small-angle scattering experiments yield the structure factor:

where $\vec{k}$ is the wave-vector of the scattered beam. A characteristic length scale $L(t)$ is usually defined as the distance at which $C(\vec{r},t)$ decays to, say, $0.2$ times its maximum value. If the domain growth is characterized by a unique length scale $L(t)$, then $C(\vec{r},t)$ and $S(\vec{k},t)$ show the dynamical scaling property Puri (2004); Puri and Wadhawan (2009): $C(\vec{r},t)=g(r/L)$; $S(\vec{k},t)=L^{d}f(kL)$. The asymptotic (large-$k$) tail of $S(\vec{k},t)$ contains information about the defects in the system. Continuous $O(n)$ spin models exhibit the generalized Porod law, with the asymptotic form: $S(k,t)\sim k^{-(d+n)}$ Porod et al. (1982); Oono and Puri (1988); Bray and Puri (1991). For $n=1$, the defects are interfaces, and the corresponding scattering function exhibits the usual Porod law: $S(k,t)\sim k^{-(d+1)}$. For $n>1$, the different topological defects are vortices ($n=2,\ d=2$), strings ($n=2,\ d=3$), and monopoles or hedgehogs ($n=3,\ d=3$). So in $d=3$, $S(k,t)\sim k^{-5}$ or $\sim k^{-6}$ depending on whether strings or monopoles dominate in the defect dynamics.

In LC mesophases such as the Sm phase, an appropriate measure of the orientational and translational order is provided by the longitudinal pair correlation function $g_{\parallel}(r_{\parallel})$ (parallel to the long axis of the LC molecules) and the transverse pair correlation function $g_{\perp}(r_{\perp})$ (perpendicular to the long axis of the LC molecules) Richter and Gruhn (2006). Evaluation of $g_{\parallel}(r_{\parallel})$ employs a cylindrical volume to probe the LC molecules aligned in the ee configuration and is given by Richter and Gruhn (2006); Cañeda-Guzmán et al. (2014):

The Heaviside step function $\theta(x)=1$ when $x\geq 0$ and $\theta(x)=0$ otherwise, $\langle\cdot\rangle$ indicates an ensemble averaging over different initial (independent) conditions, $h$ is the cylinder height used to discretize the volume, $r_{ij,\parallel_{i}}=|\bm{r}_{ij,\parallel_{i}}|=|\bm{r}_{ij}\cdot\bm{u}_{i}|$ is the center-of-mass separation along the director of molecule $i$ (the director for molecule $j$ of the pair could also be considered in this $\delta$-function evaluation as the value remains almost the same in the ee configuration). $r_{ij,\perp_{i}}=|\bm{r}_{ij,\perp_{i}}|=|\bm{r}_{ij}-\bm{r}_{ij,\parallel_{i}}|$ and $r_{ij,\perp_{j}}=|\bm{r}_{ij,\perp_{j}}|=|\bm{r}_{ij}-\bm{r}_{ij,\parallel_{j}}|$ are the corresponding transverse separations from $\bm{u}_{i}$ and $\bm{u}_{j}$. The quantity $g_{\parallel}(r_{\parallel})$ probes the average orientation of the LC molecules and the layering or smecticity in the system.

Evaluation of $g_{\perp}(r_{\perp})$ employs hollow, concentric cylinders to probe LC molecules aligned in the ss configuration and is given by Richter and Gruhn (2006); Cañeda-Guzmán et al. (2014):

In the above equation, $\delta L_{\perp}$ represents the thickness of the hollow cylinder and as for the $g_{\parallel}(r_{\parallel})$ case, we can use the director for molecule $j$ of the pair in the $\delta$-function evaluation since the value remains almost the same in the ss configuration. The quantity $g_{\perp}(r_{\perp})$ probes the translational structure about the LC molecules and their arrangement within layers.

Let us first discuss the kinetics of domain growth following a quench from the I phase to a temperature $T^{\ast}=2.5$ in the Nm phase. The left panel of Fig. 3, shows representative snapshots from the time evolution of the configurational structure for (a) $\kappa^{\prime}=5.0$, $t=8192$; (b) $\kappa^{\prime}=5.0$, $t=53248$; and (c) $\kappa^{\prime}=10.0$, $t=53248$. For clarity, we have shown only a $20^{3}$ corner of the entire box. These corners on average consist of about $N=2400$ particles. The right panel shows the corresponding top surface ($d=2$ cross-sections). There is emergence and growth of orientational Nm order, and the energy anisotropy does not affect the coarsening phenomenon. We characterize the morphologies and their texture by evaluating the correlation function and the structure factor. These are obtained by a coarse-graining procedure in which the system is divided into non-overlapping sub-boxes of size ($3.0\sigma_{0}$)³. The sub-box size is carefully chosen to ensure that each one contains about 8 to 10 particles. The continuum LC configurations are thus mapped onto a simple cubic lattice of size ($32\sigma_{0}$)³. The relevant order parameter is the orientational order parameter defined in Eq. (8), i.e., $\overline{P_{2}(\cos\theta_{i})}$. Here $\theta_{i}$ is the angle made by a molecule located in the $i$-th box with the global director $\bm{n}$ of the system (determined by taking average of orientations for all the $N$ particles present in the system), and the overline implies an average over all the particles in the sub-box.

Fig. 4(a) shows the scaled correlation function, $C(r,t)$ vs. $r/L(t)$ at $t$ = 8192, 16384 and 32768 for $\kappa^{\prime}$ = 5.0 and at $t=32768$ for $\kappa^{\prime}$ = 2.5 and 10.0. The data exhibits dynamical scaling as well as super-universality with respect to the energy anisotropy $\kappa^{\prime}$. The dynamical scaling property demonstrates that the coarsening patterns are statistically self-similar in time. The property of super-universality indicates that the morphologies in the Nm phase are independent of the energy anisotropy. This is expected because the GB energy between a pair of ellipsoids in the ss arrangement depends only on $\kappa$ and $\nu$. A log-log plot of the corresponding scaled structure factors, $S(k,t)L(t)^{-3}$ vs. $kL(t)$ is shown in Fig. 4(b). In the asymptotic large-$k$ limit, the structure factor follows the generalized Porod law: $S(k,t)\sim k^{-5}$, indicating scattering off string defects Bray et al. (1993); Birdi et al. (2020). Next, we study the time-dependence of the domain size. Fig. 4(c) shows the variation of $L(t)$ vs. $t$ on a log-log scale for $\kappa^{\prime}$ = 2.5, 5.0 and 10.0 respectively. It can be clearly observed that after an initial transient, the evolving systems are consistent with the $t^{1/2}$ growth regime (LAC law) characteristic of systems with non-conserved order parameter. The system size and time scales of our simulation are sufficient to establish the LAC domain growth law for the Nm phase of the GB model although there is onset of finite-size effects at late times.

We now come to the primary focus of our paper: kinetics of domain growth in the Sm mesophase. The coarsening is initiated by a quench from $T^{\ast}=6.0$ (I phase) to $T^{\ast}=1.0$ (Sm phase). Fig. 5 shows the prototypical evolution morphologies for: (a) $\kappa^{\prime}=5.0$, $t=4096$; (b) $\kappa^{\prime}=5.0$, $t=98304$; and (c) $\kappa^{\prime}=10.0$, $t=98304$. As in Fig. 3, we have shown only a $20^{3}$ corner of the entire box, consisting of about 2400 particles on average. The frames on the right show the corresponding top surface. It is interesting to note the initial onset of nematicity at the earlier time ($t=4096$) followed by the development of smecticity at later time ($t=98304$). Additionally, as observed in Fig. 5(c), the smecticity is significantly enhanced for $\kappa^{\prime}=10.0$.

To characterize the Sm order, Fig. 6 shows for the specified values of $t$ and $\kappa^{\prime}$: (a) the longitudinal pair correlation function $g_{\parallel}(r_{\parallel}/\sigma_{0})$ vs. $r_{\parallel}/\sigma_{0}$ evaluated using Eq. (11) and (b) the transversal pair correlation function $g_{\perp}(r_{\perp}/\sigma_{0})$ vs. $r_{\perp}/\sigma_{0}$ evaluated using Eq. (12). Notice that in Fig. 6(a), the early-time behaviour at $t=4096$ predominantly exhibits a single peak at $3\sigma_{0}$ - the length of the LC molecule. It is characteristic of Nm order with molecular alignment along an average direction or the director. As time evolves, there is emergence of newer peaks at approximately $6\sigma_{0}$, $9\sigma_{0}$, etc. as the system coarsens. This signifies development of long-range longitudinal order or layers with an inter-layer spacing of $\sim 3\sigma_{0}$. (The average separation between the peaks can be an approximate measure of the inter-layer spacing $d_{l}$. In Fig. 6(a), $d_{l}\approx 2.5\sigma_{0}$ due to interdigitation of the neighbouring layers.) Notice that increase in $\kappa^{\prime}$ reduces the intensity variation between the first and second peaks implying enhancement of smecticity. (Recall that perfect translational order is characterized by peaks of equal intensity in the pair distribution function). On the other hand, $g_{\perp}(r_{\perp}/\sigma_{0})$ vs. $r_{\perp}/\sigma_{0}$ in Fig. 6(b) exhibits peaks around multiples of $\sigma_{0}$ - the width of the ellipsoidal LC molecule. The splitting of the peak at $r_{\perp}\simeq 2.0\sigma_{0}$, a signature of BOO within the layers, is clearly seen at later times. Further, the translational order within the layers is not long-ranged as consecutive peaks have decreasing intensity. These characteristics, along with the non-zero value of the hexatic order parameter $\langle\mathcal{C}_{6}\rangle$ in Fig. 2(b), confirm the presence of the SmB-H phase. The intra-layer BOO is not affected by $\kappa^{\prime}$.

The observations from Fig. 6 suggest the following scenario for coarsening of the SmB-H phase: it is a two-timescale process, with the onset of Nm order followed by SmB-H order. To confirm this, we have also performed a similar study for the Nm$\rightarrow$SmB-H quenches and evaluated these distribution functions. Except for the emergence of multiple peaks at earlier time in the longitudinal pair correlation function (it’s now a one-timescale process as Nm order is already present and hence with time, the system exhibits only layering with BOO), the behavior of the correlation functions is qualitatively similar to that observed in Fig. 6 and hence we do not present them separately. Quantitatively, the peaks in both the functions have much higher intensity compared to the corresponding peaks in Fig. 6, indicating faster layering as well as development of BOO. We do not present these data sets to prevent repetition. The two-time-scale process is the most significant outcome from our study of the GB model. It is further reiterated by the growth laws which will be discussed shortly.

We next focus on the scaling functions that describe the time-dependent morphologies. To evaluate them, we follow the same coarse-graining approach, but now the system is divided into non-overlapping sub-boxes of size ($6.0\sigma_{0}$)³, which maps the system onto a simple cubic lattice of size ($16.0\sigma_{0}$)³. This procedure gives us a continuous order parameter field and eliminates any molecular-level anisotropies. We have changed the system size to ensure that even at later times, the sub-boxes contain about 60 to 70 particles. In Fig. 7(a), we plot $C(r,t)$ vs. $r/L(t)$ for three specified values of $t$ and $\kappa^{\prime}$. The typical scalar nematic order parameter $\overline{P_{2}(\cos\theta_{i},t)}$ is determined at each site $i$ of this discretized lattice (as for the nematic regime earlier) and subsequently the standard probes are evaluated. The data sets for $\kappa^{\prime}=5.0$ at different times neatly collapse onto a single master function, showing that the scaling regime has been reached. This data collapse indicates the existence of dynamical scaling. Furthermore, the excellent data collapse for different $\kappa^{\prime}$ values at time $t=98304$ suggests that the scaling functions are robust with respect to the anisotropy in energy. Fig. 7(d) shows the corresponding scaled structure factor, $S(k,t)L(t)^{-3}$ vs. $kL(t)$ for a range of $t$ and $\kappa^{\prime}$ values. There is an excellent data collapse, confirming both dynamical scaling and super-universality. The tail decays as $S(k,t)\sim k^{-5}$ due to the presence of string defects.

It is also possible to determine the correlation between the Sm layers by evaluating the translational correlation function $C_{l}(r,t)$. This can be done by evaluating the translational order parameter, ${\lvert\overline{\tau(\bm{r_{i}},t)}\rvert}={\lvert\overline{\exp(2i\pi\bm{r_{i}}\cdot\bm{n}/d_{l})}\rvert}$. Here, $\bm{r_{i}}$ is the position vector of the LC molecule (located in the sub-box with index $i$), $\bm{n}$ is the global director of the system and the layer spacing $d_{l}$ is evaluated from the average separation between the peaks in the longitudinal pair correlation function at a given $t$. The overline implies an average over all the particles located in the sub-box with index $i$. Fig. 7(b) shows the scaled translational correlation function $C_{l}(r,t)$ vs. $r/L_{l}(t)$ for the specified values of $t$ and $\kappa^{\prime}$. The corresponding structure factors $S_{l}(k,t)L_{l}(t)^{-3}$ vs. $kL_{l}(t)$ are shown in Fig. 7(e). The high quality of the data collapse again confirms dynamical scaling and super-universality of the scaling functions. The structure factor tail exhibits the Porod decay, $S_{l}(k,t)\sim k^{-4}$, characteristic of scattering off sharp interfaces arising between the Sm layers.

We have also determined the correlation function $C_{H}(r,t)$ using the bond-orientational order parameter ${\lvert\overline{\psi_{6}(\bm{r_{i}},t)}\rvert}$. Here, $\bm{r_{i}}$ is the position vector of the LC molecule located in the sub-box with index $i$, and the overline implies an average over all the particles located in the sub-box with index $i$. Fig. 7(c) shows the scaled correlation function $C_{H}(r,t)$ vs. $r/L_{H}(t)$ for the specified values of $t$ and $\kappa^{\prime}$. The corresponding structure factors $S_{H}(k,t)L_{H}(t)^{-3}$ vs. $kL_{H}(t)$ are shown in Fig. 7(f). The good data collapse re-confirms dynamical scaling and super-universality of the scaling functions.

Finally, we study the domain growth laws for the SmB-H mesophase. Fig. 8(a) shows $L(t)$ vs. $t$, Fig. 8(b) shows $L_{l}(t)$ vs. $t$ and Fig. 8(c) shows $L_{H}(t)$ vs. $t$ on a log-log scale for the three typical values of $\kappa^{\prime}$. (There are small differences in the pre-factors of the data sets, but they are not evident on the log-log scale.) The data sets exhibit an initial LAC growth regime $\sim t^{1/2}$ characteristic of the Nm phase and then a crossover to a slower $t^{1/4}$ growth regime. The dashed lines have been shown for reference. (Larger system sizes and longer simulation times will be required to remove the finite size effects observed at late times.) These observations emphasize the two-time-scale scenario identified in Fig. 6, and are the second novel aspect of our study.