Leandro B.
Krott†
a,
Thiago
Puccinelli
b and
José Rafael
Bordin†
*b
aCentro de Ciências, Tecnologias e Saúde, Campus Araranguá, Universidade Federal de Santa Catarina. Rua Pedro João Pereira, 150, CEP 88905120, Araranguá, SC, Brazil. E-mail: leandro.krott@ufsc.br
bDepartamento de Física, Instituto de Física e Matemática, Universidade Federal de Pelotas. Caixa Postal 354, CEP 96001-970, Pelotas, RS, Brazil. E-mail: thiagoponogueira@gmail.com; jrbordin@ufpel.edu.br
First published on 6th May 2024
Geometrical constraints offer a promising strategy for assembling colloidal crystal structures that are not typically observed in bulk or under 2D conditions. Core-softened colloids, in particular, have emerged as versatile chemical building blocks with applications across various scientific and technological areas. In this study, we investigate the behavior of a core-softened model confined between two parallel walls. Employing molecular dynamics simulations, we analyze the system's response under extreme confinement, where only one or two layers of colloids are permitted. The system comprises particles modeled by a ramp-like potential confined within slit nanoslits created by two flat, purely repulsive walls with a lateral side L separated by a distance Lz. Through a systematic analysis of the phase behavior as Lz increases, or as the system undergoes decompression, for different values of L, we identified a mono-to-bilayer transition associated with changes in the colloidal structure. In the monolayer regime, we observed solid phases at lower densities than those observed in the 2D case. Importantly, we demonstrated that confinement at specific Lz values, allowing particle arrangement into two layers, can lead to the emergence of the square phase, which was not observed under monolayer or 2D conditions. By correlating thermodynamic, translational, and orientational ordering, as well as the dynamics of this confined colloidal system, our findings offer valuable insights into the utilization of geometrical constraints to induce and manipulate structural changes.
For example, polymeric core–shell particles typically comprise a solid, often inorganic core such as silica or polystyrene, enveloped by a soft polymer shell like polyethylene glycol (PEG) or poly(N-isopropylacrylamide) (PNIPAM). Similarly, microgel particles typically feature a soft, deformable shell encasing a crosslinked polymer hard core. The softness of these colloids imparts unique properties due to the simple interplay between the hard core and soft shell length scales in their interactions. Moreover, the straightforward control over the characteristic length scales makes them excellent chemical building blocks for obtaining a wide variety of bottom-up self-assembled morphologies in nanostructured materials.5–10
In this regard, it is noteworthy how confinement within narrow slits significantly alters the behavior of colloidal systems, introducing unique dynamics and properties.11–15 In these constrained spaces, colloidal and solute particles experience restricted motion and altered interactions,11,16,17 leading to distinctive phenomena compared to bulk systems.18 Confinement can induce structural changes, such as the formation of ordered and disordered layers of particles.14,19 Moreover, the presence of surfaces in close proximity can influence particle–surface interactions, resulting in adsorption or depletion effects. By adjusting the slit size, shape, and surface properties, it is possible to manipulate phase transitions and enhance colloidal stability,12,17,20–23 impacting the rheological and mechanical properties of the system.24–26 Therefore, by precisely tuning the dimensions and surface properties of the confining environment, researchers can exert control over colloidal assembly, manipulate nucleation kinetics, and control the formation of crystal polymorphs, enabling the design of novel materials with tailored functionalities for various applications.17,27–32
Many core-softened (CS) approaches have been proposed to model these complex fluids, including ramp-like, square and ramp shoulders, exponential, Hertzian, multi-Hertzian, and Gaussian potentials.33–42 For instance, in the work by Grillo and co-authors,5 a Hertzian potential was employed to reproduce the distinct patterns observed in soft colloids at liquid–liquid interfaces. Interestingly, they found that a key ingredient to observe a large variety of micropatterns is to first immobilize two layers of colloids onto a solid substrate one after the other. This approach allowed them to obtain structures that would not be observed if the particles were assembled in a single step at the water–oil interface. Another colloidal system with competitions between two conformations is composed of pure or grafted PEG colloids, polymer-grafted nanoparticles, and star polymers, and many studies have indicated that these systems are well described by ramp-like potentials.43–48 This approach assumes that the repulsive interactions between these macromolecules arise from a gradual increase in repulsion as these colloids approach each other, often mimicking the behavior of a ramp. This model accounts for the steric repulsion between polymer chains, which dominates the interactions at short distances, while also incorporating the effects of the particle's hard core.49 Originally proposed by Hemmer and Stell50 to study polyamorphism, many ramp-like potentials, continuous or discontinuous, have been employed to study systems with two characteristic length scales and their unique behavior, which includes solid and liquid polymorphism, polyamorphism, and water-like anomalies.33,51–65 Thus, ramp-like models are able to reproduce what is observed in experiments, such as colloidal clustering, formation of self-assembled mesophases, layering packing, and even synthesis of new structures.66–72
De Oliveira and colleagues have demonstrated that a continuous ramp-like shape can be constructed by combining a Lennard-Jones potential with a Gaussian well. This simple and purely repulsive model was initially used to reproduce and study water-like anomalies in bulk73,74 and under confinement.75,76 Recently, we revisited this 3D bulk model,77 uncovering the solid phases for this model and demonstrating that an order-disordered transition in the solid phase is related to the water-like anomalies. In the 2D limit, Cardoso et al.78 have shown that this model can reproduce complex tessellations observed in experiments. For instance, in the work by Rey et al.,79 it was observed that poly(N-isopropylacrylamide) nanogels can assemble into distinct patterns at liquid–liquid interfaces, including two triangular (or hexagonal) phases, one densely packed and one less dense, as well as stripe and square phases. The ramp-like model employed in this work also exhibits two distinct triangular phases, which we refer to as low and high-density triangular phases, respectively, and a stripe pattern phase,78 as observed in the aforementioned experimental investigation.79 Triangular and stripe patterns were also experimentally observed in the self-assembly of superparamagnetic spheres.66 Additionally, the polymer-like fluid, where colloids form linear aggregates, obtained in experiments of poly(ethylene glycol) chains grafted onto poly(styrene) particles, was reproduced by the potential. Also, this spherically symmetric approach led to a Kagome lattice phase, similar to the one reported experimentally for patchy colloids.80,81 Additionally, the ramp-like model led to reentrant melting phases between the distinct crystals, which are related to a second and a third region of water-like anomalies, and later shown that these extra anomalous regions vanish as the system transitions from a 2D-like, layered regime under confinement to the 3D limit.19 This, along with inspiration from experimental works,5,66,79 led to another question: if confinement can induce anomalous regions, can it induce new crystalline morphologies in the system?
To answer this question, we carried out molecular dynamics simulations of this CS fluid confinement in narrow slits – going from the limit where only one (2D limit) layer of colloids can fit and isothermically expanding to wider slits where two layers can fit. Our findings show how it is possible to obtain specific morphologies by tuning the wall separations, including a square crystal that was not observed in the 2D case. Combining thermodynamic, dynamic, and structural analysis, we demonstrate how the distinct structures can be obtained by simply controlling the nanoslits size. This paper is organized as follows: in Section 2, we present the model, the methods, and simulation details; in Section 3, we show and discuss the results; and in Section 4, we present the conclusions.
(1) |
Using the parameters u0 = 5ε, c02 = 1.0, and r0/σ = 0.773,74,77,83 the potential 1 exhibits the ramp-like shape shown in Fig. 1 (solid purple line). As de Oliveira and co-authors have shown,84 both the real and the imaginary branch of the instantaneous normal mode (INM) spectra of this potential have a pronounced bimodality which must be connected with two different length scales – unlike simple liquids, as Lennard Jones fluids, who have only one scale. The first length-scale, correspondent to the hard-core, is located near r1 = 1.2σ, where the force has a local minimum,84 while the longer length scale, the soft corona, is located at r = 2.0σ, where the fraction of imaginary modes of the INM spectra has a local minimum and a maximum is expected in the radial distribution function84 – a schematic depiction of a polymer grafted nanoparticle as a CS particle can be seen in the inset of Fig. 1. Also, this specific set of parameters allowed us the observe distinct patterns in the 2D system, as open and close triangular lattices, stripes, kagome lattice,78 which are also obtained experimentally.5,66,80 The cutoff radius for the interaction is rc = 3.5σ.
In our system, N = 1000 particles were confined between parallel, smooth walls, whose interaction potential is purely repulsive, given by R6 potential, as follow19,23,85
(2) |
The simulations were conducted using a custom-made program in the canonical ensemble, which maintains a fixed number of particles N, volume V, and temperature T. A quasi-2D geometry was achieved by placing fixed confining walls with an area of L × L in the xy-plane separated by a distance Lz. Subsequently, systems with distinct separations, ranging between 2.0 and 4.4, were simulated. In these scenarios, the systems formed either one or two layers, referred to as monolayers or bilayers, respectively.
The particle density ρ is calculated as
(3) |
The temperature was held constant using the Nosé–Hoover thermostat with a coupling parameter Q = 2, and the equations of motion were integrated using the velocity Verlet algorithm. Simulations were conducted with 1 × 106 time steps for equilibration, followed by an additional 1 × 106 time steps to obtain averages and results. A time step of δt = 0.001 was considered. Equilibrium was monitored by analyzing pressure and energy as functions of time. Our simulations were performed along the isotherm T = 0.075, where both the 2D and 3D cases are in the solid phase.77,78
The structure was analyzed by means of the lateral radial distribution function (LRDF) g‖(r‖), defined as
(4) |
(5) |
The cumulative two-body entropy also is a property obtained from the LRDF and a useful tool to analyze the particles longer-range structure. It is defined as90
(6) |
The orientational order for particles inside a slab of thickness δz = σ was characterized with the bond orientational order parameter Ψl,
(7) |
(8) |
Once the strong confinement in the z-direction prevents diffusion along this axis, we proceeded to evaluate the dynamics within the xy-plane using the lateral diffusion coefficient,
(9) |
Furthermore, due to the slab geometry, it is important to compute the perpendicular pressure, P⊥, defined by
All quantities used in our simulations were given in LJ units.82 Some examples of reduced units are: reduced distance r* = r/σ, reduced temperature T* = kbT/ε, reduced pressure P* = σ2P/ε and reduced time . For simplicity, the * will be omitted.
The larger plates exhibit a smooth behavior, as highlighted by the black triangles connected by the solid lines for L = 90.0, with the pressure decreasing as Lz increases. This behavior is reinforced by its derivative, shown in Fig. 2(b) depicted by the black triangles. On the other hand, at intermediate values of L, the perpendicular pressure initially decreases with Lz expansion down to a minimum, then increases to a maximum and decreases once again. This loop-like behavior is highlighted in Fig. 2(a) by the blue curves, with the curve for L = 60.0 being highlighted with blue diamonds connected by a solid line.
Although this behavior contradicts intuition since one would anticipate a decrease in P⊥ as the nanoslit size grows, it is striking how this pattern resembles a van der Waals loop observed in phase transitions. At lower L's, not only one but two loops were observed, as evidenced by the red square curves in Fig. 2(a), with the case L = 40.0 being highlighted with a solid line. The derivative of this curve is shown in Fig. 2(b). Similar van der Waals loops were observed in a previous work for the present model related to structural transitions between two and three layers.94 Therefore, we can expect that here it can be also related to it. For this reasons, a useful quantity to analyze at this point is the density profile along the confined direction.
Starting with the case where L = 90.0, we can see from Fig. 3(a) that the system assumes a single-layer arrangement, or monolayer, when it is strongly compressed. This is indicated by the blue curves and the inset. As the separation, Lz, increases, the layer thickness also increases, as shown by the dashed black curves. For the wider nanoslits, the density profile assumes a bulk-like profile, represented by the red curves. Although for the larger separation, Lz = 4.4, there are two layers, the fact that ρ(z) ≠ 0 in the space between them indicates that the particles can change from one layer near a wall to another layer, exhibiting fluid-like behavior.
We can also check the colloidal particles arrangement by looking at the LRDF for this wall size, Fig. 3(d): essentially the same fluid-like g‖(r‖) is observed for all values of separation. The long-range translational ordering isn't affected, as indicated by the convergence of the cumulative two-body entropy in Fig. 3(g). Then, although there is a change from a 2D-like layering to a bulk-like layering, there is no change in the lateral structure, but as a consequence, the P⊥ × Lz curve is smooth. Also, the density range for this plate separation, L = 90, goes from ρ = 0.124 when Lz = 2.0 to ρ = 0.036 for the wider plate separation, Lz = 4.4. This means that the system is in a low density regime.
Moving on to the case, where the plate size L = 60.0 corresponds to a intermediate density regime, with 0.277 ≥ ρ ≥ 0.081 along the decompression curve. The Fig. 3(b) indicates that the arrangement goes from quasi-2D under high compression to a thicker monolayer as Lz increases – similar to the L = 90.0 case. However, the scenario changes when the separation is large enough to accommodate two layers. The black dashed curves show the separations where the layer splits into two well-defined layers – or the mono-to-bilayer transition. Interestingly, this transition from a monolayer to a bilayer coincides with the loop in P⊥ illustrated in Fig. 2(a).
Moreover, the LRDF indicates that the loop for L = 60.0 is also associated with a change in the colloids' structure. As depicted in Fig. 3(e), when the particles are arranged in a monolayer fashion, they display a fluid-like g‖(r‖) with short-range ordering, transitioning to a solid-like g‖(r‖) with longer-range ordering for wider nanoslits. This behavior might seem counter-intuitive, since one would expect a loss in ordering under decompression. However, this can be attributed to the properties of core-softened colloids.77,78 We illustrate in Fig. 3(e) that, as Lz expands, the occupancy in the first length scales diminishes while growing in the second length scale. The larger slit allows the particles to conform to a well-defined structure dominated by this greater length scale, while higher compression – or density – frustrates the system, leading to a less ordered structure dominated by the closer length scale. As indicated by the |C(r‖)| in Fig. 3(h), the solid-like phase exhibits long-range ordering, spanning to a distance R = L/2, whereas the fluid-like phase displays short-range ordering, converging at short distances.
Finally, the L = 40.0 case,where the system is in a higher density regime, 0.625≥ρ≥0.183, is depicted in Fig. 3(c), also exhibits a loop corresponding to the transition from a mono to a bi-layer arrangement. On the other hand, the second loop is associated with an increase in the peak in ρ(z) as Lz expands, as observed by comparing this case with the L = 60.0 case. Additionally, under the condition of small Lz, the colloids adopt a solid-like structure with long-range ordering, as shown in Fig. 3(f). Upon decompression, the structure transits to a fluid-like ordering, as indicated by the dashed red lines – this order–disorder transition corresponds to the first loop in Fig. 2(a). The second loop emerges when the fluid structure changes, as evident from the comparison of the dashed and solid curves in Fig. 3(f), resulting in a loss of longer-range ordering. This is supported by the cumulative two-body entropy. In the monolayer region, Fig. 3(i) suggests long-range ordering. Conversely, in the bilayer regime, two distinct behaviors are observed: both indicate short-range ordering, but one (the solid lines, corresponding to the wider slits) is shorter. This transition between two disordered phases corresponds to the second loop in Fig. 2.
To relate the fluid structure and mobility, we can juxtapose the behavior of the translational and orientational order parameters τ and Ψl with the self-diffusion coefficient D, illustrated in Fig. 4(a), (b), and (c), respectively. For instances with two loops, such as L = 40.0, we observe that τ diminishes with the plates' separation – as expected, given by the decreasing density. However, we can delineate three regimes: the ordered one, within the monolayer structure, and two exhibiting lower ordering in the bilayer region. Notably, the orientational pair ordering, Ψ2, also escalates. Evaluating the values of Ψl in the 2D-limit (Lz = 2.0), they resemble those observed in the Kagome phase for the 2D scenario,78,95 and, as depicted in Fig. 5(a), the system indeed adopts a Kagome monolayer crystal structure. Consequently, upon decompression, it transitions into two distinct clusterized fluids: at intermediate Lz values, worm-like clusters emerge (the initial fluid phase in the bilayer regime, Fig. 5(b)), and as Lz increases, these clusters break into smaller, dimer-like formations, as shown in Fig. 5(c) – for clarity, particles in different layers are color-coded (blue or red). Additionally, the diffusion coefficient, illustrated in Fig. 4(c), indicates that in the Kagome phase, particles are practically immobile – akin to a solid phase. Subsequently, they exhibit slow diffusion in the worm-like fluid phase and heightened mobility in the dimer-like fluid conformation. Consequently, we are able to illustrate the scenario along the P⊥ × Lz curve: the initial loop corresponds to a transition from a Kagome crystal to a worm-like fluid, coinciding with the colloids transitioning from a monolayer to bilayer arrangement. The subsequent loop denotes a transition between two clusterized fluids: worm-like and dimer-like.
Fig. 4 (a) Translational (τ) and (b) orientational (Ψ), respectively, ordering parameters and (c) lateral diffusion coefficient as function of Lz. |
Once the scenario for this curve was revealed, we applied this procedure to all the other curves in diagram 2(a). Along the curves with one loop, which are associated with the mono-to-bilayer transition, we observed distinct morphologies evolving in those transitions. For instance, let's compare the L = 50.0 and L = 60.0 cases. They correspond to densities ranges 0.400 ≥ ρ ≥ 0.118 and 0.277 ≥ ρ ≥ 0.081, respectively. The parameter τ, displayed in Fig. 4(a), indicates that along both curves, the system shifts from a less ordered to a more ordered morphology. Likewise, the diffusion coefficient D, illustrated in Fig. 4(c), signals a transition from higher to lower mobility. However, for L = 50.0, the less ordered phase indeed exhibits a value of τ that indicates high ordering, whereas for L = 60.0, τ aligns more closely with structured fluids.74,77,96 Additionally, we can observe in Fig. 4(b) that for L = 50.0, the system boasts a high value of Ψ2 in the monolayer regime, corresponding to a Stripe phase, as depicted in Fig. 5(d). As the slit expands allowing the formation of two layers, Ψ6 rises and dominates the particle orientation, as expected for a triangular phase. In our previous works,78,95 we referred to this as a Low-Density Triangular (LDT) phase, given that the particle separation equals the soft shell scale, as seen in Fig. 3(e). The snapshot in Fig. 5(d) showcases the LDT lattice in both layers – here, we connected the blue colloids separated by distances smaller than rb = 3.0 for better visualization of the LDT structure. Therefore, for L = 50.0, the loop in Fig. 2(a) constitutes a transition between two ordered phases: a Stripe phase in the monolayer regime to an LDT phase in the bilayer structure.
When it comes to the mono-to-bilayer transition for L = 60.0, the scenario is different. Here, there is a monolayer fluid for narrower slits, characterized by lower translational and orientational ordering and higher diffusion, and an ordered phase inside wider slits with D ≈ 0, as indicated by Fig. 4(a), (b), and (c) respectively. Even more surprisingly, the higher value of Ψ4 in the bilayer ordered phase indicates a square phase that was not observed in the 2D case.78 In fact, the snapshots in Fig. 5(f) and (g) support the fluid-to-square crystal transition associated with the mono-to-bilayer transition. Furthermore, it is worth noting that the square phase was only observed when the system was arranged in two layers – no single-layer case exhibited this crystalline structure. This observation, coupled with the absence of this phase in bulk 2D, reinforces that this morphology is induced by confinement in our model. These findings align with those obtained by Grillo and colleagues:5 the frustration induced by the confinement between two walls led to the square morphology that was not observed in the 2D case.78 Also, is relevant to remember that a square phase was not observed in the simulations for this model, but can be found in another CS models.97
Finally, even for the case where there is no loop in the P⊥ × Lz curve, we can observe effects of the mono-to-bilayer transition. For the case where L = 90.0, when the system is in the fluid phase, it is observed that the translational and orientational order parameters are smaller compared to the ordered phases, with respect to the value of Lz. Also, we would expect a decrease in fluid ordering and an increase in diffusion as the density decreases – or, in our case, as the slit widens. However, we can see in Fig. 4(a), (b), and (c) that for L = 90.0, this transition corresponds to an increase in τ and Ψ2 and a decrease in D. Although the model employed in this work, depicted by eqn (1), does exhibit water-like anomalies – indicating an increase in ordering and decrease in diffusion under expansion – we believe that confinement also plays a role, particularly in fluid dynamics. Fig. 4(c) displays a sharp decrease in D when there is enough space to start the formation of the second layer – as shown in Fig. 5(h) and (i), the snapshots for the system at Lz = 3.2 and Lz = 3.6, where the system transitions from higher to lower diffusion regimes – and where the density distribution width is changing from the 2D limit to the bilayer regime, as depicted in Fig. 3(a).
When using the entire procedure throughout this work, such as analyzing the density histogram, RDF, cumulative two-body entropy, orientational and translational order parameters, diffusion coefficient and snapshots applied to all values of L simulated, or all P⊥ × Lz, we were able to construct the ρ × LZ state diagram shown in Fig. 6.
While it serves as a qualitative state diagram, Fig. 6 illustrates how confinement in narrow spaces impacts the behavior of soft colloids. Within the monolayer region, depicted with a blue background, the system displays the same phases as those observed in the 2D case: initially a fluid (abbreviated to Fl in the diagram), followed by the LDT phase as the density increases. Subsequently, a clusterized fluid phase separates the ordered LDT and Stripes phases. In the high-density regime, a Stripes and Kagome coexistence region, denoted as StK in the phase diagram, emerges, and eventually, the Kagome structure is achieved. An intriguing observation within the monolayer region is that as Lz increases, the confined ordered phases can be observed at lower densities compared to the 2D case. For example, the Kagome crystal was observed at densities ρ ≈ 0.40 for Lz = 2.6, lower than in the 2D limit,78 and similar observations were made for the Stripes and LDT phases.
In between the monolayer and bilayer regimes lies the mono-to-bilayer fluid region, delineated by the magenta background. When arranged in two layers, represented by the red region in Fig. 6, we observe the LDT and square phases at low densities. Both structures transit to a Stripe pattern at higher densities, but while the square and Stripes phase is separated by a cluster phase, the LDT and Stripes coexist in the TtS region. Finally, akin to the monolayer regime, the Stripes phase transitions to the Kagome phase, with the StK coexistence region demarcating the phases.
Our findings revealed distinct features in the P⊥ × Lz curves, particularly evident during the transition from one to two layers. We explored these transitions in detail for intermediate plate sizes, such as L = 50.0 and L = 60.0, where we observed transitions from stripes to low-density triangular (LDT) and from fluid to solid square phases, respectively. Additionally, for smaller plate sizes like L = 40.0, we observed a transition from Kagome to worm-like fluid. Some curves, such as L = 40.0, exhibited a second loop, indicating a transition from a worm-like fluid to one characterized by smaller, dimer-like clusters. Conversely, for larger plates such as L = 90.0, no evidence of a structural transition was observed. However, a notable decrease in the diffusion coefficient D was observed as the wall separation increased, corresponding to the emergence of the second layer and indicating a significant impact of the layering transition on the system's dynamics.
To demonstrate the utility of confinement in controlling the aggregation process for this model, we constructed the ρ × Lz diagram. In the monolayer regime, our findings indicated that open crystal phases, including LDT, Stripes, and Kagome, could be obtained at lower densities with increasing slit width. In the bilayer regime, a notable emergence of a square ordered phase, not observed in the 2D case, was detected. Those finding highlights how the extreme confinement inside narrow slits can be employed to stabilize quasi-2D structures and even obtain new arrangements.
Footnote |
† These authors contributed equally to this work. |
This journal is © The Royal Society of Chemistry 2024 |