Short, partially complementary, single-stranded (ss)DNA strands can form nanostructures with a wide variety of shapes and mechanical properties. It is well known that semiflexible, linear dsDNA can undergo an isotropic to nematic (IN) phase transition and that sufficiently bent structures can form a biaxial nematic phase. Here, we use numerical simulations to explore how the phase behavior of linear DNA constructs changes as we tune the mechanical properties of the constituent DNA by changing the nucleotide sequence. The IN-phase transition can be suppressed in so-called DNA “nunchakus”: structures consisting of two rigid dsDNA arms, separated by a sufficiently flexible spacer. In this paper, we use simulations to explore what phase behavior to expect for different linear DNA constructs. To this end, we first performed numerical simulations exploring the structural properties of a number of different DNA oligonucleotides using the oxDNA package. We then used the structural information generated in the oxDNA simulations to construct more coarse-grained models of the rod-like, bent-core, and nunchaku DNA. These coarse-grained models were used to explore the phase behavior of suspensions of the various DNA constructs. The approach explored in this paper makes it possible to “design” the phase behavior of DNA constructs by a suitable choice of the constituent nucleotide sequence.

In his seminal 1949 paper,1 Onsager predicted that thin, rod-like colloids could undergo a purely entropy-driven transition from an isotropic liquid to a nematic liquid crystal (LC).

Entropy-driven phase transitions of suspensions of non-spherical, hard colloidal particles have been the subject of many experimental, theoretical and numerical studies; see, e.g., Refs. 2–5. In addition, many experimental and simulation studies focused on thermotropic molecular liquid crystal formers amongst which “bent-core” or “banana-shaped” molecules garnered particular interest because they display a rich phase behavior that includes different smectic phases. Such liquid crystals are of interest in view of their potential applications in non-linear optics and display technology.6–11 It is challenging to synthesize colloidal bent-core mesogens with sufficiently narrow size- and angle-distributions such that a clear liquid crystal (LC) phase transition could be observed.12 An important advance was made by Fernańdez-Rico et al.,13 who succeeded in synthesizing banana-shaped colloids with a low polydispersity for which the bending angle could be controlled continuously. Fernańdez-Rico et al.13 performed confocal-microscopy studies on suspensions of their banana-shaped colloids to explore the various LC phases directly as a function of the bending angle. Their observations could be reproduced in variational mean-field theory and Monte Carlo simulations.10 

An attractive approach to design perfectly monodisperse, bent-core systems with no polarity or attractive inter-particle potentials is the use of relatively short, single-stranded DNA chains. Such chains had earlier been used to create self-assembling DNA nanostars14–16 and a wide range of DNA-origami structures.17–22 The large difference in the persistence length of double- and single-stranded DNA has also been used to explore linear, double-stranded (ds)DNA as liquid crystal formers.23–27 However, in experiments, a key challenge for using DNA lies in the design of thermodynamically stable DNA sequences that can then organize into the desired LC-phase. A practical problem is that designed DNA sequences tend to be relatively short (for cost reasons) and that for the typical aspect ratios of such short DNA constructs high concentrations are needed to reach the IN-phase transition. As a consequence, it is costly to prepare many different mesogenic DNA constructs.

Salamonczyk et al.25 introduced a tailor-made, coarse-grained model to specifically emulate the smectic phase transition of DNA nunchakus they studied in experiments. There they used a 20 thymine long spacer to connect two rigid dsDNA duplexes, which were represented by hard cylinders connected by two beads via appropriate springs. This choice allowed the two rigid arms to align parallel to each other, thus promoting smectic alignment. However, this model did not reflect the dependence of the bending angle between the two arms as a function of the ssDNA linker-length. Here, we show how to use oxDNA 2.0,28,29 a semi-coarse-grained simulation model that provides structural and thermodynamic properties of any DNA construct, to design nearly rigid rods and corresponding rigid, bent-core, and flexible “nunchakus” mesogens and then test these in the laboratory. Subsequently, we transposed our computationally expensive oxDNA model into a further coarse-grained bead-spring model developed by Xing et al.16 for DNA nanostars using Large-scale Atomic/Molecular Massively Parallel Simulator (LAMMPS).30 This approach allowed simulating larger systems of DNA mesogens and provides a better predictor of the experimental phase behavior of such systems. We identified individual phases by measuring the diffusion coefficients of the mesogens along the different coordinates. While the literature exists on the dynamic behavior of mesogens within liquid crystal phases,31,32 this method was not used previously to identify the nature of the smectic phase.

We base our coarse-grained bead-spring model16 of rigid rods and nunchakus with different flexibilities on the molecular structure and interaction parameters obtained from the semi-empirical, more detailed oxDNA calculations28,33 of DNA sequences, which we tested experimentally. The specific sequence and detailed experimental and oxDNA modeling are given in the supplementary material.

The nunchaku molecules are formed of two rigid rods connected by a flexible linker, as depicted in Fig. 1. Each rod is modeled as a series of beads that have a diameter of 2 nm, corresponding to ∼6 base pairs of dsDNA.34 These beads are connected by harmonic springs and kept linear by an angular potential with a minimum at 180°. All molecules considered here are shorter than the persistence length of double-stranded DNA, which is ∼50 nm.35 Therefore, the approximation of near-perfect rigidity (and linear rods) is valid and justified further in the supplementary material.

FIG. 1.

Both schematic representations of our coarse-grained bead-spring model and the oxDNA calculation below the melting temperature of the DNA sequences used.

FIG. 1.

Both schematic representations of our coarse-grained bead-spring model and the oxDNA calculation below the melting temperature of the DNA sequences used.

Close modal

Two separate cases for the central angular potential connecting the two rods were considered, allowing for bent-core and flexible mesogens. Bent-core molecules represent the limiting case of a completely rigid connector with the minimum of the angular potential set to θ = 150°. This allows for small fluctuations around θ, thus giving all mesogens approximately the same opening angle. Such molecules would require a carefully designed linking sequence tailored to the opening angle desired, and we propose that techniques from DNA nanotechnology, such as those pioneered by Siavashpouri et al.,36 would be best placed to realize this experimentally.

In the flexible case, the energy scale of the angular potential is reduced while the minimum is kept at 180° so that the opening angle of the mesogens can vary considerably. This choice was based on our oxDNA analysis of the flexibility provided by single-stranded DNA that links the two arms in the nunchakus. These oxDNA simulations, detailed in the supplementary material and stated also in previous work,16 suggest that an ssDNA linker of minimum length of 6 bases (roughly one bead size) provides a wide probability distribution of θ as shown in Fig. 1. In the same figure, we also show the angle-probability distribution for the same sequence with no ssDNA flexibility–this quasi-hard rod system was used to calibrate our simulations with regard to isotropic to nematic phase transitions as a function of aspect ratio.

The model was implemented in LAMMPS30 using reduced units based on the fundamental units of mass m, the Lennard-Jones (LJ) energy scale ϵLJ, the LJ length scale σLJ, and the Boltzmann constant kB. Unless otherwise stated, these units are used in all figures throughout. These also define a characteristic time scale τ=(mσLJ2/ϵLJ)1/2. The LJ quantities are defined by the form of the inter-molecular interaction potential for which a shifted, cut-off potential (potted in Fig. 1) is used,
Uij=4ϵLJσLJrij12σLJrij6+ϵLJ,rij<rc=21/6σLJ,
(1)
where Uij = 0 otherwise. In this Weeks–Chandler–Andersen (WCA) style potential, a cut-off distance rc = 21/6σLJ was chosen to represent the purely repulsive excluded volume interactions between the beads, ensuring that all phase transitions observed are truly entropy-driven. Neighboring beads are connected by the harmonic potential,
Vbond=Kbond(rr0)2,
(2)
where r0 is the equilibrium bond distance and Kbond is the stiffness of the harmonic bond. We set r0 to 0.96σLJ and Kbond to 500ϵLJ/σLJ2 throughout, allowing for minimal disturbances around the equilibrium distance. As the beads themselves have a radius 0.5rc (equivalent to 0.56σLJ), this gives a small overlap between neighboring beads. The angular potential is given by
Vangle=Kangle(θθ0)2,
(3)
where θ0 is the equilibrium bond angle and Kangle is the spring constant that controls the allowed distribution of the angle. We set θ0rod=180° within each rod and θ0bend=150° for the bent-core molecules. We chose Kangle = 500ϵLJ/rad2 throughout, except for the central angle in the flexible nunchakus, where it is given as 0.1ϵLJ/rad2 (a non-zero value is used when defining this bond’s properties to ensure numerical stability in the simulation). This corresponds to the ssDNA, represented by a separate bead in the center of the molecule with differing mechanical properties.

All molecular dynamics (MD) simulations presented here were performed in LAMMPS.30 Simulations were conducted on a system of 1000 particles with a time step of 0.005τ, unless otherwise stated. Furthermore, the simulations were performed within an oblong box defined by the Cartesian axes under periodic boundary conditions.37 The aspect ratio of this box may be varied to support phase formation for anisotropic mesogens and evaluate correlation functions over longer lengths without increasing the size of the system.

The system was initially configured in a dilute, isotropic state, prepared by random placement and orientation of the mesogens in a simulation box using a Monte Carlo algorithm. Any molecules that did not fit in the simulation box or those that overlapped with the existing units were discarded, and the placement procedure was repeated until all particles were added. When studying high volume fractions, simulations were initiated from a perfectly ordered square crystalline phase with all molecules aligned along a common axis. Care was taken to ensure stability of this ordered phase by confirming molecules did not overlap and that the internal energy was conserved after equilibration.

An isenthalpic ensemble was used to expand/contract the size of the simulation box, hence varying the volume fraction of the mesogen system, while a microcanonical (N, V, E) ensemble was used to thermally equilibrate the system at each volume fraction studied. For the isenthalpic ensemble, time integration was evaluated using the Nosé–Hoover barostat,38,39 natively implemented in LAMMPS40 with a damping time τ.

Simulation runs consisted of multiple stages alternating between these two ensembles to sweep through different system densities: contraction stages varied the volume fraction while equilibration stages were used for data collection. Typically 2 × 104 steps (each of duration 0.005τ) were simulated in the isenthalpic ensemble, followed by 2 × 106 steps in the microcanonical ensemble, to allow the system to reach equilibrium at this volume fraction. Static system properties to characterize the equilibrium phase of the system were sampled after each equilibration period. To ensure the stability of the system during these periods, energy and temperature conservation (at T = 0.5ϵLJ/kB) were verified over a range of timescales.

Our linear, rod-like DNA mesogens display traditional liquid crystal phase transitions with characteristic nematic and smectic order parameters. The nematic order parameter Sn is given by
Sn=P2(cosφ)=32cos2φ12,
(4)
where P2 denotes the second Legendre polynomial41 and φ is the angle between the molecular axis and the system director. This is non-trivial to calculate when the system director is not known (i.e., in the absence of an external field). Therefore, we used the approach taken by Frenkel and Eppenga42–further details are provided in the supplementary material.
Similarly, the smectic order parameter is defined as
Ss=1Nj=1Nexp2πLirjn̂
(5)
for N molecules in layers of periodicity L perpendicular to the nematic director n̂; rj denotes the center-of-mass position of the jth molecule.43 The value of L that maximizes Ss corresponds to the smectic layer spacing or the characteristic length scale of positional order along the axis n̂.
It is also instructive to introduce a length-dependent orientational order parameter to verify that orientational-ordering of molecules is truly long-ranged, occurring across all length scales, and is not simply due to short-range steric effects. We consider an lth rank, pair-wise correlation function, where gl(r) gives the correlation between the orientation of two particles (determined by the angle between their molecular axes) separated by a distance r,
gl(r)=Pl(ûiûj)δ(rijr)δ(rijr),
(6)
where ûi is the director for molecule i and rij is the separation between a given pair of molecules i and j (determined from their molecular center-of-mass positions).44 In the disordered phase, this pair-correlation function decays to zero, while in the ordered phase it decays to the square of the orientational order parameter,45 
limrgl(r)=Pl2.
(7)

It is worth noting here that the maximum separation between particles is half the size of the simulation region due to the periodic boundary conditions46 and that this will vary over the duration of the simulation. To determine order over longer length scales without increasing the volume of the system, we varied the aspect ratio of the simulation region and sampled correlations along the long axis of the box.

The explicit calculation of gl(r) scales as O(N2) and requires storing the angle between every pair of molecules, limiting the resolution possible.47 We may avoid this problem by expanding the function in terms of spherical harmonics and instead sum (2l + 1) contributions from the spherical harmonics for each molecule. Further details of this approach are provided in the supplementary material.

A system of rigid rods was used to benchmark the analysis methods in comparison to theoretical predictions derived by Onsager1 and, in particular, the isotropic–nematic phase transition, which has been verified computationally through both Monte Carlo48,49 and molecular dynamics simulations50,51 as well as experimentally.52–54 

For rigid rods with an aspect ratio L/D = 10, we confined the possible critical volume fraction for this lyotropic transition to the range 0.39 < ϕ < 0.44, as observed in Fig. 2, achieving good agreement with Onsager’s prediction of ϕ = 0.4. This analysis was also repeated for longer rods with an aspect ratio of 16 with simulations identifying the critical volume fraction in the range 0.23 < ϕ < 0.26 in good agreement with the predicted value of ϕ = 0.25.

FIG. 2.

The evolution of the volume fraction and the nematic order parameter over the timescale of the simulation for a system of 1000 rigid rods with aspect ratio 10. The discrete change in Sn (red line) occurs at about ϕ ∼ 0.4, denoting the transition to the nematic phase, with the system depicted explicitly on either side of this transition. Note that the contraction steps, where the volume fraction is increased, are not of equal duration and so do not correspond to equal changes in the system volume; rather, they are chosen to highlight the phase transition.

FIG. 2.

The evolution of the volume fraction and the nematic order parameter over the timescale of the simulation for a system of 1000 rigid rods with aspect ratio 10. The discrete change in Sn (red line) occurs at about ϕ ∼ 0.4, denoting the transition to the nematic phase, with the system depicted explicitly on either side of this transition. Note that the contraction steps, where the volume fraction is increased, are not of equal duration and so do not correspond to equal changes in the system volume; rather, they are chosen to highlight the phase transition.

Close modal

We also consider the phase behavior upon expansion from a perfectly ordered, crystalline state. This has two advantages: it allows access to higher volume fractions that are not easily accessible through the shrinking procedure and provides verification of the phase transitions previously observed. Ensuring that a novel phase is in true equilibrium has long been a challenge for liquid crystal simulators; however, non-equilibrium effects will manifest themselves in a hysteresis of the phase transition (variation in the critical volume fraction dependent on the direction of the transition) and can, therefore, be identified through this method.

The isotropic phase formation was observed in the region 0.38 < ϕ < 0.41 in good agreement with the previous simulations, confirming the equilibrium nature of this phase transition. The higher volume fractions accessed at the start of the simulation also suggested the existence of a smectic phase with a sharp change in the smectic order parameter around ϕ = 0.6; however, this has not been investigated further due to the wealth of the literature on this already.55,56

These techniques were subsequently applied to the nunchaku molecules introduced in Sec. II A. We consider first the flexible linker case, which we have verified provides an accurate coarse-grained model of the DNA-nunchakus using oxDNA.29 The simulation molecules were lengthened to include 15 beads instead of 10 as minimal order was observed at the lower rod aspect ratios.

The linker flexibility was sufficient to ensure a quasi-isotropic initial angle distribution in the range 60° < θ < 180° with smaller angles excluded due to steric repulsion between the rods. Unlike previous studies,25 this prohibits the nonphysical possibility of the “fully bent” (θ = 0°) conformation where both rods are parallel in agreement with the prior oxDNA results. The angle distribution we observed from the detailed oxDNA model indicated the low possibility of the fully bent situation. More information from the oxDNA measurement can be found in the supplementary material.

Despite this flexibility, a strong preference for the linear configuration (θ = 180°) was observed with a quasi-nematic phase (Sn = 0.6) formed at volume fractions ϕ > 0.4. However, no uniform, system-wide director was observed with the preferential orientation varying across the simulation region shown in Fig. 3. It is possible that periodic variation of the director occurs on a length scale greater than the size of the simulation region, but the computational resources required to simulate larger systems exceeded those available, preventing further characterization of this phase.

FIG. 3.

The kernel density estimate of the opening angle distribution as the volume fraction is reduced for completely flexible nunchakus. Plotted for a system of 103 particles with the distribution sampled every 1.5 × 106 time steps. Note the formation of a preferential angle at late times, corresponding to the formation of an ordered phase at high volume fractions, the final form of which is depicted in the inset.

FIG. 3.

The kernel density estimate of the opening angle distribution as the volume fraction is reduced for completely flexible nunchakus. Plotted for a system of 103 particles with the distribution sampled every 1.5 × 106 time steps. Note the formation of a preferential angle at late times, corresponding to the formation of an ordered phase at high volume fractions, the final form of which is depicted in the inset.

Close modal

To force the formation of more exotic LC phases with anisotropic (non-linear) mesogens, we also considered bent-core rods with a fixed opening angle. We considered opening angles θ = 60°, 90°, 120°, 150°, 165° and found that larger opening angles favored nematic-like phase formation (with a maximal order parameter Sn = 0.62), but smaller angles did not form obvious alternate ordered phases in agreement with previous surveys of this range of opening angles.11,57 In this section, we focus on θ = 150°. These findings are consistent with the previous literature on the subject–early work in the area adopts a bend angle of θ = 140°,59–60 while more recent studies have tended to favor values closer to θ = 150°.62–63 

Using a pair-wise orientational correlation function (the calculation of which is discussed further in the supplementary material), we confirm that the order observed in this fixed angle system is indeed long-range with a non-decaying correlation at long length scales.

This method is limited by the use of a single vector along the molecular axis to define the orientation of the molecule; two vectors are required to uniquely specify the orientation of a single nunchaku molecule. Using the molecular axis alone only accounts for quasi-nematic order and does not consider any biaxial or twisted nematic substructure. We therefore also considered additional orientational correlation functions for the bisector of each molecule and normal vector to the plane of the nunchaku. The first Legendre polynomial is used to detect changes in sign of the bisector direction, as would be expected in a herringbone structure.

However, no statistically significant periodic components in the correlation function were observed over the length scale of the simulation region. A visual inspection of the inset in Fig. 4 suggests that the length of the simulation region is approximately half the full period of the repeating “twist” in the ordered phase, suggesting the possible existence of a twisted nematic phase in this system. However, further research on larger and more computationally expensive simulation regions with greater length scales than the periodic twist is required to verify these phases.

FIG. 4.

Orientational correlation function over time (as the volume fraction is reduced) for a system of 1000 nunchaku molecules with a fixed opening angle of 150 × 102 deg and sampled every 7 × 105 time steps. Note the formation of an ordered phase at high volume fractions (late times) with sustained long-range order (i.e., no decay in the orientational correlation function). The maximum value of particle separation is determined by the size of the simulation box and so shrinks over time.

FIG. 4.

Orientational correlation function over time (as the volume fraction is reduced) for a system of 1000 nunchaku molecules with a fixed opening angle of 150 × 102 deg and sampled every 7 × 105 time steps. Note the formation of an ordered phase at high volume fractions (late times) with sustained long-range order (i.e., no decay in the orientational correlation function). The maximum value of particle separation is determined by the size of the simulation box and so shrinks over time.

Close modal
While the phase behavior of rigid rods is well studied, much less is known about the dynamic properties of the individual mesogens within these quasi-ordered phases. Nevertheless, dynamic studies through MD simulations have enjoyed recent popularity in both computational and experimental work.65–66 Here, we study the dynamic properties of the nunchaku system and demonstrate the application of dynamic properties to static phase identification. In particular, we focus on the diffusion coefficient D, defined by the “power-law” equation for mean-squared displacement (MSD) from an initial center-of-mass position x0 over time t in n dimensions,67,
x(t)x02=2nDαtα,
(8)
which reduces to the pure-diffusive case when the power exponent α = 1. Otherwise, the process is characterized as subdiffusive (α < 1) or superdiffusive (α > 1).68 

Displacements were only sampled over equilibration (constant volume) stages of our simulation, and the effect of the periodic boundary conditions was explicitly accounted for by offsetting additional displacements from boundary crossings. Linear regression analysis was used to validate the value of α for each run over a variety of volume fractions in the dilute limit, where each particle may exist in a non-overlapping free volume of rotation (a sphere circumscribed around the molecule). For rigid rods with an aspect ratio of 10, this occurs at ϕd = 0.015. We found an average power of α = 0.97 ± 0.03 for the nunchaku molecules and α = 1.00 ± 0.04 for the rigid rods, as expected for diffusive behavior.

As the volume fraction was increased, nunchakus display highly subdiffusive behavior with a reduction in the value of α to 0.12 ± 0.02 at ϕ = 0.69. We also observe anisotropy in the coordinate-wise displacement in non-isotropic, ordered phases, as demonstrated in Fig. 5, which considers three different phase regimes for a system of rigid rods. Here, the Cartesian axes are used as simulations are configured in a crystalline phase with molecules orientated along the y axis and smectic layers forming in the xz plane. Phase transitions for this system were previously identified in the regions 0.58 < ϕ < 0.62 for the smectic–nematic and 0.38 < ϕ < 0.43 for the nematic–isotropic transition.

FIG. 5.

Root mean-squared displacements against time in a system of 1000 rigid rods with an aspect ratio 10 for a range of volume fractions ϕ. The leftmost pair of plots correspond to the smectic phase with restricted motion between xz layers along the y axis. The center pair give the nematic phase with preferential displacement along the y axis. Finally, the rightmost pair give the isotropic phase with isotropic diffusion and no preferential direction. Note that no significant differences appear between the x- and z-directions in any phase as these directions are equivalent in the phase structure.

FIG. 5.

Root mean-squared displacements against time in a system of 1000 rigid rods with an aspect ratio 10 for a range of volume fractions ϕ. The leftmost pair of plots correspond to the smectic phase with restricted motion between xz layers along the y axis. The center pair give the nematic phase with preferential displacement along the y axis. Finally, the rightmost pair give the isotropic phase with isotropic diffusion and no preferential direction. Note that no significant differences appear between the x- and z-directions in any phase as these directions are equivalent in the phase structure.

Close modal

At the highest volume fractions (ϕ ≥ 0.62), the diffusion along the y axis is significantly reduced due to the restricted motion between layers in the smectic phase. Below this, a nematic phase is observed with increased diffusivity along the molecule’s director, corresponding to an increased intercept in logarithmic space, as observed for the y axis in the central panels of Fig. 5. The degree of dynamic anisotropy here, quantified by the ratio of diffusion coefficients Dy/Dx, is 2.12 for ϕ = 0.58. At the lowest volume fractions (ϕ ≤ 0.38), there is no preferred direction in the system and motion is isotropic and truly diffusive (as α = 1).

Phase formation may be directly observed through deviations in these coordinate-specific diffusion coefficients. Within the flexible nunchaku system, isotropic diffusion is observed at low volume fractions in Fig. 6, but the y-coordinate diffusion coefficient is severely reduced beyond this with a maximal anisotropy of Dy/Dx = 0.19 at ϕ = 0.58. This strongly suggests the presence of a smectic phase below ϕ = 0.45 in agreement with the prediction of 0.40 < ϕ < 0.44 obtained from the smectic order parameter.

FIG. 6.

Coordinate diffusion coefficients for a microcanonical ensemble of flexible nunchaku at a range of volume fractions ϕ. The low volume fraction region on the left of the graph corresponds to the isotropic phase with no variation between coordinate axes. In contrast, the smectic phase in the high volume fraction gives rise to anisotropy in the coordinate diffusion coefficients with reduced diffusivity perpendicular to the smectic layers in the y axis direction.

FIG. 6.

Coordinate diffusion coefficients for a microcanonical ensemble of flexible nunchaku at a range of volume fractions ϕ. The low volume fraction region on the left of the graph corresponds to the isotropic phase with no variation between coordinate axes. In contrast, the smectic phase in the high volume fraction gives rise to anisotropy in the coordinate diffusion coefficients with reduced diffusivity perpendicular to the smectic layers in the y axis direction.

Close modal

Dynamic properties, such as directional diffusion coefficients, may be utilized to detect novel phases that are otherwise difficult to characterize. Building on theoretical predictions by Camp et al.,69,70 we determined that an opening angle θ = 120° would be optimal for biaxial phase formation. We were able to demonstrate the existence of a biaxial smectic phase through the measurements of directional diffusion coefficients, as depicted in Fig. 7. The system was initially configured as a crystalline lattice of linear molecules before the minima in the opening angle were shifted to form the desired opening angle. The initial orientation of the bisector axis is random, and so, initially, no preferential bisector direction was observed. Subsequent equilibration (2 × 106 steps) allowed relaxation into the smectic biaxial state, which was then observed through repeated simulations (2 × 106 steps) to transition back into the traditional smectic and isotropic phases upon expansion.

FIG. 7.

Root mean-squared displacements against time in a system of 1000 bent-core mesogens with an opening angle 120° over a range of volume fractions ϕ. As previously, restricted diffusion is observed along the director (previously the y axis) between smectic layers. Anisotropy is now observed within these layers at the highest volume fraction ϕ = 0.76, suggesting a biaxial smectic phase with a common aligned bisector. Diffusion is reduced in the direction of the normal vector compared to the bisector vector for these mesogens.

FIG. 7.

Root mean-squared displacements against time in a system of 1000 bent-core mesogens with an opening angle 120° over a range of volume fractions ϕ. As previously, restricted diffusion is observed along the director (previously the y axis) between smectic layers. Anisotropy is now observed within these layers at the highest volume fraction ϕ = 0.76, suggesting a biaxial smectic phase with a common aligned bisector. Diffusion is reduced in the direction of the normal vector compared to the bisector vector for these mesogens.

Close modal

This demonstrates the existence of a biaxial smectic phase for this system of DNA mesogens and the power of dynamic property simulation in novel phase identification.

We considered the phase behavior of DNA nunchaku molecules, consisting of two sections of dsDNA connected by a short section of ssDNA. We introduced a coarse-grained model of DNA nunchaku particles, formed of two rigid rods connected by a (semi-) flexible linker, with soft-core, purely repulsive potential interactions. We also introduced a simpler “rigid rod” model, without the flexible linker, to verify the phase identification techniques used here.

Our simulations of rigid rods were consistent with Onsager’s prediction of the existence of a first-order “entropic” phase transition between the isotropic and nematic phases of slender hard rods across a wide range of aspect ratios. In comparison to Onsager’s predictions of a transition volume fraction ϕ = 0.4 for rods with aspect ratio L/D = 10, we were able to confine the measured critical volume fraction within the range 0.39 < ϕ < 0.44. This value along with the equilibrium nature of this transition was verified through simulations prepared in an initial crystalline phase and subject to a stepwise expansion of the simulation box, eliminating the possibility of phase hysteresis.

The application of the same techniques to the nunchaku system gave evidence for the formation of an ordered, quasi-nematic phase in this system in both the semiflexible (Sn = 0.58) and bent-core (Sn = 0.62) configurations. The fixed rigidity configuration, a more realistic representation of the DNA system, demonstrated clear evidence for an entropy-driven phase transition through the reduction in configurational entropy associated with the formation of a preferred angle. The pair-wise orientational correlation function was also used to validate the presence of truly long-range order.

Finally, we demonstrate that the measurement of dynamic properties provides a little-studied, alternative method for identifying phase transitions and characteristic symmetries of such systems in comparison with phases previously identified in both the rigid-rod and nunchaku systems. In particular, we considered the formation of anisotropic phases through the variation in coordinate diffusion coefficients, such as in the smectic phase where the diffusion coefficient within the ordered layers is up to a factor of 5.15 greater than that between them. We extend this characterization of diffusive properties beyond the arbitrary Cartesian coordinates imposed by the simulation region to allow the identification of novel phases where the global director or symmetry class is not a priori known and must be inferred from the symmetry properties of the system. Using this approach to directional diffusion within the nunchaku system, we also demonstrate the formation of a novel biaxial phase not previously observed in DNA-based systems. We suggest that this method can be used in further identification of novel liquid crystal phases with directional anisotropy and is particularly applicable in previously uncharacterized systems where the symmetry properties are unknown, or there is not an appropriate order parameter tailored to the transition of interest.

Through a greater understanding of the phase behavior of DNA nanoparticles, we may obtain better insight into the design of more complex mesogen shapes and thus develop new, entropy-driven self-assembled structures.71,72 These techniques may also have applications in biophysics, photonics, structural biology, and synthetic biology.74–75 

The supplementary material details the oxDNA model, including the specific base pair sequences and subsequent coarse-graining process to generate the mesogen model used throughout. It also outlines computational processes for the nematic and positional order parameters, specifically deriving the spherical harmonic form for the latter.

K.G. acknowledges funding from the Streetly Fund for Natural Sciences at Queens’ College, Cambridge. J.Y. and R.L. acknowledge the Cambridge Trust and China Scholarship Council (CSC), and D.A.K. acknowledges the UK Engineering and Physical Sciences Research Council Ph.D. Studentship (Award No. 1948692) for financial support. E.E. acknowledges support of the Research Council of Norway through its Centres of Excellence funding scheme, Project No. RCN 262644.

The authors have no conflicts to disclose.

Kit Gallagher: Investigation (lead); Writing – original draft (equal). Jiaming Yu: Investigation (equal); Supervision (supporting); Writing – original draft (supporting). David Andrew King: Investigation (supporting). Ren Liu: Investigation (supporting). Erika Eiser: Conceptualization (lead); Supervision (lead); Writing – review & editing (equal).

The data that support the findings of this study along with all simulations and analysis codes are openly available at https://github.com/KCGallagher/LC_Project. Data for work conducted in oxDNA are available separately.76 

1.
L.
Onsager
, “
The effects of shape on the interaction of colloidal particles
,”
Ann. N. Y. Acad. Sci.
51
,
627
659
(
1949
).
2.
V. J.
Anderson
and
H. N. W.
Lekkerkerker
, “
Insights into phase transition kinetics from colloid science
,”
Nature
416
,
811
815
(
2002
).
3.
R.
Eppenga
and
D.
Frenkel
, “
Monte Carlo study of the isotropic and nematic phases of infinitely thin hard platelets
,”
Mol. Phys.
52
,
1303
1334
(
1984
).
4.
M.
Engel
,
P. F.
Damasceno
,
C. L.
Phillips
, and
S. C.
Glotzer
, “
Computational self-assembly of a one-component icosahedral quasicrystal
,”
Nat. Mater.
14
,
109
116
(
2015
).
5.
M.
Dijkstra
, “
Entropy-driven phase transitions in colloids: From spheres to anisotropic particles
,” in
Advances in Chemical Physics
, edited by
S.
Rice
and
A.
Dinner
(
John Wiley & Sons, Inc.
,
2015
), Vol.
156
, pp.
35
71
.
6.
H.
Takezoe
and
Y.
Takanishi
, “
Bent-core liquid crystals: Their mysterious and attractive world
,”
Jpn. J. Appl. Phys.
45
,
597
625
(
2006
).
7.
J.
Etxebarria
and
M.
Blanca Ros
, “
Bent-core liquid crystals in the route to functional materials
,”
J. Mater. Chem.
18
,
2919
(
2008
).
8.
R. A.
Reddy
and
C.
Tschierske
, “
Bent-core liquid crystals: Polar order, superstructural chirality and spontaneous desymmetrisation in soft matter systems
,”
J. Mater. Chem.
16
,
907
961
(
2006
).
9.
M.
Chiappini
,
T.
Drwenski
,
R.
van Roij
, and
M.
Dijkstra
, “
Biaxial, twist-bend, and splay-bend nematic phases of banana-shaped particles revealed by lifting the ‘smectic blanket
,’”
Phys. Rev. Lett.
123
,
068001
(
2019
).
10.
M.
Chiappini
and
M.
Dijkstra
, “
A generalized density-modulated twist-splay-bend phase of banana-shaped particles
,”
Nat. Commun.
12
,
2157
(
2021
).
11.
P.
Kubala
,
W.
Tomczyk
, and
M.
Cieśla
, “
In silico study of liquid crystalline phases formed by bent-shaped molecules with excluded-volume type interactions
,”
J. Molec. Liquids
367
,
120156
(
2021
).
12.
Y.
Yang
,
H.
Pei
,
G.
Chen
,
K. T.
Webb
,
L. J.
Martinez-Miranda
,
I. K.
Lloyd
,
Z.
Lu
,
K.
Liu
, and
Z.
Nie
, “
Phase behaviors of colloidal analogs of bent-core liquid crystals
,”
Sci. Adv.
4
,
eaas8829
(
2018
).
13.
C.
Fernández-Rico
,
M.
Chiappini
,
T.
Yanagishima
,
H.
de Sousa
,
D. G. A. L.
Aarts
,
M.
Dijkstra
, and
R. P. A.
Dullens
, “
Shaping colloidal bananas to reveal biaxial, splay-bend nematic, and smectic phases
,”
Science
369
,
950
(
2020
).
14.
N. C.
Seeman
, “
Nucleic acid junctions and lattices
,”
J. Theor. Biol.
99
,
237
247
(
1982
).
15.
S.
Biffi
,
R.
Cerbino
,
F.
Bomboi
,
E. M.
Paraboschi
,
R.
Asselta
,
F.
Sciortino
, and
T.
Bellini
, “
Phase behavior and critical activated dynamics of limited-valence DNA nanostars
,”
Proc. Natl. Acad. Sci. U. S. A.
110
,
15633
15637
(
2013
).
16.
Z.
Xing
,
C.
Ness
,
D.
Frenkel
, and
E.
Eiser
, “
Structural and linear elastic properties of DNA hydrogels by coarse-grained simulation
,”
Macromolecules
52
,
504
512
(
2019
).
17.
P. W. K.
Rothemund
, “
Folding DNA to create nanoscale shapes and patterns
,”
Nature
440
,
297
302
(
2006
).
18.
Y.
Dong
,
Z.
Yang
, and
D.
Liu
, “
DNA nanotechnology based on i-motif structures
,”
Acc. Chem. Res.
47
,
1853
1860
(
2014
).
19.
Y. J.
Cha
,
M.-J.
Gim
,
K.
Oh
, and
D. K.
Yoon
, “
Twisted-nematic-mode liquid crystal display with a DNA alignment layer
,”
J. Inf. Disp.
16
,
129
135
(
2015
).
20.
L.
Wang
,
A. M.
Urbas
, and
Q.
Li
, “
Nature-inspired emerging chiral liquid crystal nanostructures: From molecular self-assembly to DNA mesophase and nanocolloids
,”
Adv. Mater.
32
,
1801335
(
2018
).
21.
N. C.
Seeman
and
H. F.
Sleiman
, “
DNA nanotechnology
,”
Nat. Rev. Mater.
3
,
17068
(
2017
).
22.
J. W.
Canary
and
H.
Yan
, “
Nadrian C. (Ned) Seeman (1945–2021)
,”
Nat. Nanotechnol.
17
,
108
(
2022
).
23.
F.
Livolant
and
A.
Leforestier
, “
Condensed phases of DNA: Structures and phase transitions
,”
Prog. Polym. Sci.
21
,
1115
1164
(
1996
).
24.
T. P.
Fraccia
,
G. P.
Smith
,
L.
Bethge
,
G.
Zanchetta
,
G.
Nava
,
S.
Klussmann
,
N. A.
Clark
, and
T.
Bellini
, “
Liquid crystal ordering and isotropic gelation in solutions of four-base-long DNA oligomers
,”
ACS Nano
10
,
8508
8516
(
2016
).
25.
M.
Salamonczyk
,
J.
Zhang
,
G.
Portale
,
C.
Zhu
,
E.
Kentzinger
,
J. T.
Gleeson
,
A.
Jakli
,
C.
De Michele
,
J. K. G.
Dhont
,
S.
Sprunt
, and
E.
Stiakakis
, “
Smectic phase in suspensions of gapped DNA duplexes
,”
Nat. Commun.
7
,
13358
(
2016
).
26.
C.
De Michele
,
G.
Zanchetta
,
T.
Bellini
,
E.
Frezza
, and
A.
Ferrarini
, “
Hierarchical propagation of chirality through reversible polymerization: The cholesteric phase of DNA oligomers
,”
ACS Macro Lett.
5
,
208
212
(
2016
).
27.
D.
Winogradoff
,
P. Y.
Li
,
H.
Joshi
,
L.
Quednau
,
C.
Maffeo
, and
A.
Aksimentiev
, “
Chiral systems made from DNA
,”
Adv. Sci.
8
,
2003113
(
2021
).
28.
T. E.
Ouldridge
,
A. A.
Louis
, and
J. P. K.
Doye
, “
Structural, mechanical, and thermodynamic properties of a coarse-grained DNA model
,”
J. Chem. Phys.
134
,
085101
(
2011
).
29.
P.
Šulc
,
F.
Romano
,
T. E.
Ouldridge
,
L.
Rovigatti
,
J. P. K.
Doye
, and
A. A.
Louis
, “
Sequence-dependent thermodynamics of a coarse-grained DNA model
,”
J. Chem. Phys.
137
,
135101
(
2012
).
30.
A. P.
Thompson
,
H. M.
Aktulga
,
R.
Berger
,
D. S.
Bolintineanu
,
W. M.
Brown
,
P. S.
Crozier
,
P. J.
in’t Veld
,
A.
Kohlmeyer
,
S. G.
Moore
,
T. D.
Nguyen
,
R.
Shan
,
M. J.
Stevens
,
J.
Tranchida
,
C.
Trott
, and
S. J.
Plimpton
, “
LAMMPS—A flexible simulation tool for particle-based materials modeling at the atomic, meso, and continuum scales
,”
Comput. Phys. Commun.
271
,
108171
(
2022
).
31.
E.
De Miguel
,
L. F.
Rull
,
M. K.
Chalam
, and
K. E.
Gubbins
, “
Liquid crystal phase diagram of the Gay-Berne fluid
,”
Mol. Phys.
74
,
405
424
(
1991
).
32.
A. D.
Rey
, “
Theory and simulation of gas diffusion in cholesteric liquid crystal films
,”
Mol. Cryst. Liq. Cryst. Sci. Technol., Sect. A
293
,
87
109
(
1997
).
33.
J. P. K.
Doye
,
T. E.
Ouldridge
,
A. A.
Louis
,
F.
Romano
,
P.
Šulc
,
C.
Matek
,
B. E. K.
Snodin
,
L.
Rovigatti
,
J. S.
Schreck
,
R. M.
Harrison
, and
W. P. J.
Smith
, “
Coarse-graining DNA for simulations of DNA nanotechnology
,”
Phys. Chem. Chem. Phys.
15
,
20395
20414
(
2013
).
34.
S.
Arnott
and
D. W. L.
Hukins
, “
Optimised parameters for A-DNA and B-DNA
,”
Biochem. Biophys. Res. Commun.
47
,
1504
1509
(
1972
).
35.
H. G.
Garcia
,
P.
Grayson
,
L.
Han
,
M.
Inamdar
,
J.
Kondev
,
P. C.
Nelson
,
R.
Phillips
,
J.
Widom
, and
P. A.
Wiggins
, “
Biological consequences of tightly bent DNA: The other life of a macromolecular celebrity
,”
Biopolymers
85
,
115
130
(
2007
).
36.
M.
Siavashpouri
,
C. H.
Wachauf
,
M. J.
Zakhary
,
F.
Praetorius
,
H.
Dietz
, and
Z.
Dogic
, “
Molecular engineering of chiral colloidal liquid crystals using DNA origami
,”
Nat. Mater.
16
,
849
856
(
2017
).
37.
D.
Frenkel
and
B.
Smit
,
Understanding Molecular Simulation: From Algorithms to Applications
(
Academic Press
,
San Diego
,
2002
), Chap. 3, pp.
32
35
.
38.
S.
Nosé
, “
A unified formulation of the constant temperature molecular dynamics methods
,”
J. Chem. Phys.
81
,
511
519
(
1984
).
39.
W. G.
Hoover
, “
Canonical dynamics: Equilibrium phase-space distributions
,”
Phys. Rev. A
31
,
1695
1697
(
1985
).
40.
W.
Shinoda
,
M.
Shiga
, and
M.
Mikami
, “
Rapid estimation of elastic constants by molecular dynamics simulation under constant stress
,”
Phys. Rev. B
69
,
134103
(
2004
).
41.
P.
de Gennes
and
J.
Prost
,
The Physics of Liquid Crystals
,
International Series of Monogr
(
Clarendon Press
,
1993
).
42.
D.
Frenkel
and
R.
Eppenga
, “
Evidence for algebraic orientational order in a two-dimensional hard-core nematic
,”
Phys. Rev. A
31
,
1776
1787
(
1985
).
43.
S.
Dussi
,
M.
Chiappini
, and
M.
Dijkstra
, “
On the stability and finite-size effects of a columnar phase in single-component systems of hard-rod-like particles
,”
Mol. Phys.
116
,
2792
2805
(
2018
).
44.
C.
Zannoni
, in
The Molecular Physics of Liquid Crystals
, edited by
G. R.
Luckhurst
and
G. W.
Gray
(
Academic Press
,
London
,
1979
), Chap. 3, p.
51
.
45.
D.
Frenkel
and
B. M.
Mulder
, “
The hard ellipsoid-of-revolution fluid
,”
Mol. Phys.
55
,
1171
1192
(
1985
).
46.
D.
Frenkel
,
B. M.
Mulder
, and
J. P.
Mctague
, “
Phase diagram of hard ellipsoids of revolution
,”
Mol. Cryst. Liq. Cryst.
123
,
119
128
(
1985
).
47.
A. K.
Soper
, “
Determination of the orientational pair correlation function of a molecular liquid from diffraction data
,”
J. Mol. Liq.
78
,
179
200
(
1998
).
48.
D.
Frenkel
,
B. M.
Mulder
, and
J. P.
McTague
, “
Phase diagram of a system of hard ellipsoids
,”
Phys. Rev. Lett.
52
,
287
290
(
1984
).
49.
S. D.
Lee
, “
A numerical investigation of nematic ordering based on a simple hard-rod model
,”
J. Chem. Phys.
87
,
4972
4974
(
1987
).
50.
M. P.
Allen
and
D.
Frenkel
, “
Observation of dynamical precursors of the isotropic-nematic transition by computer simulation
,”
Phys. Rev. Lett.
58
,
1748
1750
(
1987
).
51.
P. J.
Camp
,
C. P.
Mason
,
M. P.
Allen
,
A. A.
Khare
, and
D. A.
Kofke
, “
The isotropic–nematic phase transition in uniaxial hard ellipsoid fluids: Coexistence data and the approach to the Onsager limit
,”
J. Chem. Phys.
105
,
2837
2849
(
1996
).
52.
K.
Kubo
and
K.
Ogino
, “
Comparison of osmotic pressures for the poly (y-benzyl-L-glutamate) solutions with the theories for a system of hard spherocylinders
,”
Mol. Cryst. Liq. Cryst.
53
,
207
227
(
1979
).
53.
R.
Oldenbourg
,
X.
Wen
,
R. B.
Meyer
, and
D. L. D.
Caspar
, “
Orientational distribution function in nematic tobacco-mosaic-virus liquid crystals measured by x-ray diffraction
,”
Phys. Rev. Lett.
61
,
1851
1854
(
1988
).
54.
S.
Fraden
,
G.
Maret
, and
D. L. D.
Caspar
, “
Angular correlations and the isotropic-nematic phase transition in suspensions of tobacco mosaic virus
,”
Phys. Rev. E
48
,
2816
2837
(
1993
).
55.
D.
Frenkel
,
H. N. W.
Lekkerkerker
, and
A.
Stroobants
, “
Thermodynamic stability of a smectic phase in a system of hard rods
,”
Nature
332
,
822
823
(
1988
).
56.
S. C.
McGrother
,
D. C.
Williamson
, and
G.
Jackson
, “
A re-examination of the phase diagram of hard spherocylinders
,”
J. Chem. Phys.
104
,
6755
6771
(
1996
).
57.
T.
Drwenski
and
R.
van Roij
, “
The effect of flexibility and bend angle on the phase diagram of hard colloidal boomerangs
,”
Mol. Phys.
116
,
2812
2822
(
2018
).
58.
R.
Memmer
, “
Liquid crystal phases of achiral banana-shaped molecules: A computer simulation study
,”
Liq. Cryst.
29
,
483
496
(
2002
).
59.
S. J.
Johnston
,
R. J.
Low
, and
M. P.
Neal
, “
Computer simulation of apolar bent-core and rodlike molecules
,”
Phys. Rev. E
65
,
051706
(
2002
).
60.
J.
Xu
,
R. L. B.
Selinger
,
J. V.
Selinger
, and
R.
Shashidhar
, “
Monte Carlo simulation of liquid-crystal alignment and chiral symmetry-breaking
,”
J. Chem. Phys.
115
,
4333
4338
(
2001
).
61.
A.
Dewar
and
P. J.
Camp
, “
Computer simulations of bent-core liquid crystals
,”
Phys. Rev. E
70
,
011704
(
2004
).
62.
A.
Vaghela
,
P. I. C.
Teixeira
, and
E. M.
Terentjev
, “
Emergence of biaxial nematic phases in solutions of semiflexible dimers
,”
Phys. Rev. E
96
,
042703
(
2017
).
63.
H. G.
Yoon
,
S.-W.
Kang
,
R. Y.
Dong
,
A.
Marini
,
K. A.
Suresh
,
M.
Srinivasarao
, and
S.
Kumar
, “
Nematic biaxiality in a bent-core material
,”
Phys. Rev. E
81
,
051706
(
2010
).
64.
F.
Gay-Balmaz
,
T. S.
Ratiu
, and
C.
Tronci
, “
Equivalent theories of liquid crystal dynamics
,”
Arch. Ration. Mech. Anal.
210
,
773
811
(
2013
).
65.
T.
Zhao
and
X.
Wang
, “
Diffusion of rigid rodlike polymer in isotropic solutions studied by dissipative particle dynamics simulation
,”
Polymer
54
,
5241
5249
(
2013
).
66.
A. D.
Rey
,
E. E.
Herrera-Valencia
, and
Y. K.
Murugesan
, “
Structure and dynamics of biological liquid crystals
,”
Liq. Cryst.
41
,
430
451
(
2013
).
67.
D.
Ernst
and
J.
Köhler
, “
Measuring a diffusion coefficient by single-particle tracking: Statistical analysis of experimental mean squared displacement curves
,”
J. Chem. Phys.
15
,
000845
849
(
2013
).
68.
R.
Metzler
and
J.
Klafter
, “
The random walk’s guide to anomalous diffusion: A fractional dynamics approach
,”
Phys. Rep.
339
,
1
77
(
2000
).
69.
P. J.
Camp
and
M. P.
Allen
, “
Phase diagram of the hard biaxial ellipsoid fluid
,”
J. Chem. Phys.
106
,
6681
6688
(
1997
).
70.
P. J.
Camp
,
M. P.
Allen
, and
A. J.
Masters
, “
Theory and computer simulation of bent-core molecules
,”
J. Chem. Phys.
111
,
9871
9881
(
1999
).
71.
E.
Barry
and
Z.
Dogic
, “
Entropy driven self-assembly of nonamphiphilic colloidal membranes
,”
Proc. Natl. Acad. Sci. U. S. A.
107
,
10348
10353
(
2010
).
72.
K.
Lin
,
J. C.
Crocker
,
V.
Prasad
,
A.
Schofield
,
D. A.
Weitz
,
T. C.
Lubensky
, and
A. G.
Yodh
, “
Entropically driven colloidal crystallization on patterned surfaces
,”
Phys. Rev. Lett.
85
,
1770
1773
(
2000
).
73.
S.
Nummelin
,
J.
Kommeri
,
M. A.
Kostiainen
, and
V.
Linko
, “
Evolution of structural DNA nanotechnology
,”
Adv. Mater.
30
,
1703721
(
2018
).
74.
F.
Praetorius
,
B.
Kick
,
K. L.
Behler
,
M. N.
Honemann
,
D.
Weuster-Botz
, and
H.
Dietz
, “
Biotechnological mass production of DNA origami
,”
Nature
552
,
84
87
(
2017
).
75.
M.
Bathe
and
P. W. K.
Rothemund
, “
DNA nanotechnology: A foundation for programmable nanoscale materials
,”
MRS Bull.
42
,
882
888
(
2017
).
76.
K.
Gallagher
,
J.
Yu
,
D. A.
King
,
R.
Liu
, and
E.
Eiser
(
2022
). “
Towards new liquid crystal phases of DNA mesogens
,”
Zenodo
. https://doi.org/10.5281/zenodo.7871147

Supplementary Material