Next Article in Journal
Factors Influencing the Rheology of Methane Foam for Gas Mobility Control in High-Temperature, Proppant-Fractured Reservoirs
Previous Article in Journal
Saturated Micellar Networks: Phase Separation and Nanoemulsification Capacity
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Energetic and Entropic Motifs in Vesicle Morphogenesis in Amphiphilic Diblock Copolymer Solutions

by
Senyuan Liu
1 and
Radhakrishna Sureshkumar
1,2,*
1
Department of Biomedical and Chemical Engineering, Syracuse University, Syracuse, NY 13244, USA
2
Department of Physics, Syracuse University, Syracuse, NY 13244, USA
*
Author to whom correspondence should be addressed.
Colloids Interfaces 2024, 8(1), 12; https://doi.org/10.3390/colloids8010012
Submission received: 15 December 2023 / Revised: 18 January 2024 / Accepted: 29 January 2024 / Published: 4 February 2024

Abstract

:
Coarse-grained molecular dynamic simulations are employed to investigate the spatiotemporal evolution of vesicles (polymersomes) via self-assembly of randomly distributed amphiphilic diblock copolymers PB-PEO (Poly(Butadiene)-b-Poly(Ethylene Oxide)) in water. The vesiculation pathway consists of several intermediate structures, such as spherical/rodlike aggregates, wormlike micelles, lamellae, and cavities. The lamella-to-vesicle transition occurs at a constant aggregation number and is accompanied by a reduction in the solvent-accessible surface area. Simulation predictions are in qualitative agreement with the mechanism of vesicle formation in which the unfavorable hydrophobic interactions between water molecules and polymer segments, along the edge of the lamella, are eliminated at the expense of gaining curvature energy. However, rod–lamella–vesicle transition is accompanied by an increase in copolymer packing density. Hence, the change in the surface area accompanying vesiculation predicted by the simulations is significantly lower than theoretical estimates. Changes in information entropy, quantified by the expectation of the logarithm of the probability distribution function of the segmental stretch parameter s, defined as the difference between the maximum and instantaneous segmental extension, are statistically insignificant along the vesiculation pathway. For rods, lamellae, and polymersomes, s follows a log normal distribution. This is explained based on the configurational dynamics of a single diblock chain in water.

1. Introduction

The self-assembly of block copolymers (BCPs) in a solution has generated renewed interest following the pioneering discoveries of vesicular morphologies, referred to as polymersomes [1,2,3,4,5,6,7,8,9,10,11,12]. While traditional ionic surfactants with hydrocarbon chains can form vesicular structures, the manipulation of their chemical environment is often required to facilitate vesiculation, e.g., with the addition of a salt or co-surfactant, which limits their biocompatibility and biomedical applications. Further, non-ionic copolymer vesicles offer superior stability due to enhanced hydrophobicity, as evidenced by their ultra-low critical micelle concentrations compared to those of traditional surfactants. Specifically, the change in chemical potential associated with the exchange of individual copolymer molecules between a molecular assembly and surrounding solvent is large enough to render kinetic time scales spanning hours to days, as observed experimentally: see chapter 20 in [13]. Further, BCPs offer robust means of size control for the vesicles via varying polymer chain compositions/lengths, the addition of secondary spacer molecules, and microfluidics-based flow-focusing techniques to produce droplets of controllable sizes: see review by Bleul et al. [1]. Polymersomes with a size that ranges between several nanometers to tens of micrometers can be routinely manufactured. They form the building blocks of complex nanostructures that offer significant promise in detergency, environmental remediation, nanomedicine, and cancer theranostics, an approach that integrates diagnosis (e.g., by imaging) and therapeutics [1,14,15,16,17], catalysis, and nanoreactor designs [18,19,20,21] as well as the development of biomaterials that mimic cellular structures and functions [22,23,24].
Experimental studies have provided a wealth of information on and insights into the phase behavior of BCPs in a solution. Further, they have offered several clues regarding the thermodynamic and kinetic pathways underlying the morphological diversity and structure transitions in such systems. However, a quantitative theoretical understanding of energetic and entropic conduits that facilitate spontaneous changes in aggregate morphology is still lacking due to many reasons. First, the real-time visualization of the self-assembly process, leading to the formation of molecular aggregates in a solution, is well beyond the spatiotemporal resolution of advanced imaging technologies. Consequently, morphology characterization is achieved often with cryogenic transmission electron microscopy (cryo-TEM) experiments, although progress in NMR imaging has made it possible to track the trajectories of individually tagged molecules in a native solution as they navigate within and across molecular aggregates [25,26]. Second, experimental observations of shape transitions are often rationalized based on elegant yet simplified geometric arguments, such as the one based on the packing parameter, which was developed to explain spherical, cylindrical, and vesicular aggregates in amphiphilic surfactants [27]. However, experiments [7,8,9,10,11,12,28,29] and a recently reported systematic coarse-grained molecular dynamics (CGMD) study by the authors [30] have shown that the self-assembly of BCPs in a solution generates extraordinarily diverse and topologically complex morphologies. The quantitative characterization of structure evolution in such systems requires a detailed understanding of the dynamics of molecule exchange between a self-assembled structure and solvent as well as the internal reorganization of molecules that is often associated with morphology transitions. In this work, such dynamics are explored in detail in the context of vesiculation in amphiphilic diblock copolymer solutions, which involves molecular aggregation in the early stages and the internal reorganization of molecules at a constant aggregation number during the latter phase of vesiculation.
A quantitative understanding of entropic (e.g., segmental stretch and the orientation of copolymer chains) and energetic (e.g., the area of hydrophobic interfaces and potential energy of solvent–monomer interactions) changes underlying morphology transitions, especially rod-to-lamella as well as lamella-to-vesicle transitions seen in experiments [28], remains inadequate. At least two pathways have been proposed to explain vesicle formation in BCP solutions. One mechanism entails solvent diffusion into spherical micelles, creating a polymer-free core region [1]. The other mechanism, which is widely accepted by scholars due to its perceived universality to phospholipids (biological membranes) and copolymers, is based on energetic arguments, i.e., a sufficiently large lamellar (bilayer) structure develops spontaneous curvature, bends, and folds into a vesicle (or polymersome). Specifically, it is envisioned that, during such a morphology change, the lamella discards the energy associated with the unfavorable hydrophobic interactions of the polymers along its edge with the surrounding solvent molecules at the expense of gaining curvature energy [1,13,27,29,31,32,33,34,35]. In experiments conducted by Fromherz et al. [35] on egg lecithin–taurochenodesoxycholate systems, electron microscopy revealed that metastable, discoid assemblies of phospholipids spontaneously closed to form vesicular morphologies. This process is thought to be initiated only when the perimeter of the lamella exceeds a threshold value or, equivalently, when the driving force provided by the unfavorable hydrophobic edge energy, which is approximately equal to the product of the area of the lateral solvent–lamella interface and interfacial tension, is sufficiently large to overcome the bending energy of the bilayer [31]. Under such circumstances, geometric arguments, e.g., conservation of surface area of the copolymer aggregate, can provide a relationship between the linear dimensions of the lamella and vesicle produced: see for example Figure 5 in Antonietti and Förster [29].
The abovementioned conceptualization of vesicle generation from bilayers is qualitatively supported by computer simulations that utilize solvent-free multi-particle collision dynamics and discrete representations of fluid membranes parameterized to manipulate particle diffusion and membrane rigidity [36,37]. Such discrete particle simulations can capture dynamic organization of fluid membranes and provide significantly more nuanced depictions of the pathways of morphology transitions compared to those offered by classical continuum mechanic theories based on energy minimization principles that govern the bending/buckling of planar sheets to form curved and closed surfaces [38,39,40,41,42]. However, they lack critical details of the organization of constituent molecules within the aggregates, which is largely dictated by solvent-mediated interactions. In the context of aqueous solutions of copolymers, distributions of segmental stretch and water within the assemblies are highly heterogeneous and greatly dependent on the chain composition, confinement (e.g., water–hydrophilic segment interface in the inner and outer layers of a unilamellar vesicle), hydrogen bonding interactions, and/or solvent mobility [30,43]. Moreover, aggregate morphologies are dynamic with irregular shapes and diffuse interfaces [30,32,33]. Such molecular-scale aspects of shape transitions are addressed in the present study.
Experiments by Chen et al. [28] show the coexistence of rodlike micelles, lamellae, and vesicles in polystyrene-polyacrylic acid solutions in dioxane/water. These authors inferred that rod-to-vesicle transition involves an intermediate lamellar state by visualizing a series of structures by cryo-TEM. The transition from a solution that contains predominantly rods to one that is dominated by vesicles was facilitated by perturbing a solution of rodlike micelles at equilibrium by increasing its water content by incremental amounts. The coexistence of rods, lamellae, and vesicles were quantified using turbidity data. Specifically, each morphology makes a distinctive contribution to the overall turbidity of the solution, which was quantified using experimental measurements. Further, a kinetic model based on a two-state transition theory in which (a). rods are flattened to form a lamella and (b). it develops curvature and closes in on itself to create a vesicle was proposed. Intermediate morphologies between spherical micelles and vesicles, such as cylindrical micelles and lamellar (discoid) and open vesicles (curved lamellar) have been observed in experiments for amphiphilic systems [44,45]. In cetyl-trimethyl-ammonium-bromide and 5-methyl salicylic acid solution, a thermally reversible transition from vesicles to wormlike micelles has been observed [46]. In phospholipid solutions containing vesicles, a spherical micelle formation can be triggered by the addition of surfactants [44,45]. The underlying pathway consists of the opening of the vesicle to form a curved bilayer followed by a bilayer-to-rodlike micelle transition and the formation of spherical micelles from rods. In our earlier molecular dynamic study [30], the co-existence of rods, bilayers, and vesicles was predicted in diblock polymer solutions. In this work, we investigate the sequence of morphology transitions leading to vesicle formation in diblock copolymer solutions using CGMD simulations.
CGMD simulations conducted by Markvoort and coworkers [47,48] on phospholipid systems predict an increase in the overall potential energy associated with a bilayer-to-vesicle transition, suggesting the existence entropic driving forces. Whether this is true in the case of a copolymer solution is unknown. Configurational entropy associated with the distributions of a chain stretch and orientation of hydrophobic/hydrophilic segments within BCP aggregates could influence their macroscopic mechanical properties. This, in turn, could have pronounced implications on shape transitions. Establishing quantitative relationships between the molecular structure and mechanical properties of BCP aggregates is beyond the scope of the present study. However, we compute probability distribution functions of a segmental stretch for rodlike, lamellar, and vesicular aggregates to glean insights into potential entropic driving forces underlying vesiculation.
The present study is built on our recently reported CGMD simulations of the spontaneous self-assembly in a model diblock copolymer (Poly(Butadiene)-b-Poly(Ethylene Oxide), abbreviated as PB-PEO) solution leading to the formation of diverse morphologies [30]. The initial conditions of our simulations correspond to a homogeneous BCP solution at a prescribed concentration, in which polymer molecules of a given chain length and composition are randomly distributed. No pre-assembly is employed, and the spatiotemporal evolution of structures is tracked in real time. In the past, mesoscopic simulations based on self-consistent field theory [49,50,51], dissipative particle dynamics [52,53,54,55,56,57,58,59,60], and CGMD simulations with pre-assembled initial conditions [32,33,61] have been employed to study the self-assembly of amphiphilic BCPs in a solution, bilayer-to-vesicle transition, and kinetics of copolymer exchange between self-assembled aggregates in BCP solutions. Further, coarse-grained molecular dynamic simulations that account for solvent-mediated hydrophobic interactions implicitly via a suitably parameterized Lennard–Jones (LJ) potential function have been used to study differences in copolymer architecture (linear vs. bottlebrush) on self-assembly in a solution [62] and the shape deformation of PEO-PS polymersomes induced by osmotic pressure stress [63]. CGMD simulations have also been used to qualitatively study the self-assembly process of PS-PAA polymers [64], phospholipids [65], and DTAB and gemini surfactants [66]. In the present study, we employ the CGMD simulation framework reported in [30], which was adapted from our previous investigations of self-assembly and shape transitions as well as shear/elongational flow dynamics and the rheology of cationic surfactant solutions [67,68,69,70,71,72,73,74]. We use the MARTINI force field to describe all bonded and non-bonded interactions within the BCP solution [75]. Our simulations capture vesiculation from an initially homogeneous copolymer solution. This process, which is facilitated through multiple intermediate stages of micellization and structure reorganization, is analyzed to provide quantitative insights into the underlying energetics and entropic motifs. The simulation details and data analytic methods are discussed in Section 2. The results are discussed in Section 3. The conclusions and outlook are offered in Section 4.

2. Materials and Methods

2.1. Methods and Software

The force field constraints of this model were based on the requirements of the MARTINI force field [75].
The stretching and bending interactions between bonded beads were modeled by weak harmonic potentials Vb and Vθ, respectively, given by Equations (1) and (2) below:
V b   ( b ) = 1 2 K b ( b b 0 ) 2
V θ   ( θ ) = 1 2 K θ ( cos ( θ ) cos ( θ 0 ) ) 2
where b0 and θ0 represented the equilibrium bond distance and equilibrium bond angle at the minimum energy configuration, respectively, and Kb and Kθ represented constants associated with bond stretching and bending energies, respectively. The non-bonded interactions between beads i and j were described by a 6–12 Lennard−Jones (LJ) potential given by
V L J r i j = C i j 12 r 12 C i j 6 r 6 .
Equation (3) can also be written as
V L J r i j = 4 ε i j σ i j r 12 σ i j r 6 ,
where
C i j ( 12 ) = 4 ε i j σ i j 12   and   C i j ( 6 ) = 4 ε i j σ i j 6 .
In the above equations, r represents the distance between beads i and j, ε i j is the depth (minimum) of the potential well (energy function), and σ i j is the distance at which V L J is zero. ε i j and σ i j are also referred to as the interaction strength and cutoff length, respectively.
MARTINI force field requires the LJ potential to be shifted smoothly to 0 between 0.9 and 1.2 nm, rather than using a hard cutoff to avoid singularities in the computation of the inter-particle forces. This is achieved by using a potential−switch function SV given by
S V r i j = C i j 12 S 12 C i j 6 S 6 ,
where
S α = 1 r α C   ;   r r s h i f t 1 r α A 3 r r s h i f t 3 B 4 r r s h i f t 4 C   ;   r s h i f t < r r c u t 0   ;   r > r c u t
α = 6   or   12 .
In Equation (7), the parameters A, B, and C are defined by
A = α α + 4 r c u t α + 1 r s h i f t r c u t α + 2 r c u t r s h i f t 2 ,
B = α ( α + 3 ) r c u t ( α + 1 ) r s h i f t r c u t α + 2 ( r c u t r s h i f t ) 3 ,
and
C = 1 r c u t α A 3 r c u t r s h i f t 3 B 4 r c u t r s h i f t 4 ,
where values of r s h i f t = 0.9 nm and r c u t = 1.2 nm are used.
The force field parameter values used in this study were identical to those reported in Liu and Sureshkumar [30], in which morphological diversity in diblock BCP solutions mimicking PEO-PB architectures of varying chain lengths, chain compositions, and concentrations was explored. They are also listed in SI as Tables S1–S3.
In the present simulations, 400,000 water beads, 40,000 antifreeze water beads [76], and 9000 polymer chains consisting of six PEO and six PB beads (14.4 wt%) were randomly placed in an initial cubical box of linear dimension 40 nm using PACKMOL 20.14.0 software [77]. The simulations were performed using the Gromacs 2020.2 MD software. The reference pressure and temperature were 1 bar and 300 K, respectively. The production runs were NPT simulations. They spanned 400 ns, which was approximately 10% longer than the time required for a stable vesicular structure to evolve from an initially random configuration of copolymer chains. More details about the simulation are provided in SI as Section S1. The molecular visualization was performed using VMD 1.9.3. Data analysis and graphing were performed using MatlabTM 2023, Microsoft 365TM, Visual C++ 2022, and OriginlabTM 2017.

2.2. Data Analytics

Discrete spatiotemporal data generated from the CGMD simulations were analyzed to identify copolymer clusters. Subsequently, their geometric features, such as monomer number N, volume V, surface area A, and monomer density ρN/V, were computed. We also calculated solvent-accessible surface area (SASA), which may be interpreted as the surface area of the water−aggregate interface. Specifically, SASA is measured as the area circumscribed by the motion of the center of a spherical probe that is slid along the periphery of the beads that constitute a cluster [78]. The algorithms utilized for the above calculations are described in detail in our earlier work [30].

2.2.1. End-to-End Distance (Q) and Body Gaussian Curvature (KG)

To quantitatively describe shape changes of the self-assembled structures, we calculated the end-to-end distance Q and body Gaussian curvature KG of the entire assembly with respect to its center of mass. Note that KG is a wholistic (body) curvature measure of the geometry of the aggregate rather than the integral of the pointwise Gaussian curvatures along its surface. This concept is illustrated in Figure 1 for two representative surfaces, one that is relatively flat along (Figure 1a–c) and the other curved (Figure 1d–f). The Cartesian (x, y, z) coordinate system used had the center of mass of the aggregate as its origin. The shortest linear dimension was along the y-axis. Figure 1b shows the collections of discrete points (“point cloud”) that lie within −2.5 nm ≤ z ≤ 2.5 nm of the structure whose mesh representation is given in Figure 1a. This set of points was projected onto the x-y plane and binned into several intervals: see Figure 1b,e. We used 20 and 60 intervals for flat and curved surfaces, respectively, as they provide sufficient discretization to calculate numerically converged results for the curvature. A structure is considered flat if its center of mass lies within its point cloud and curved if otherwise. The curvature in the x-y plane is calculated as the inverse of the radius of the circle (Figure 1c) that fits the centers of mass of the beads in each bin (red dots in Figure 1b) by using a robust fitting algorithm developed by Pratt [79]. The data for the first two (four) and last two (four) bins for flat (curved) structures were discarded to minimize the uncertainty in the calculations arising from large fluctuations of the edges. The curvature was considered positive (negative) if the center of the best fit circle lies above (below) the center of mass of the structure. The end-to-end distance of the point cloud Qp was calculated as the magnitude of the vector connecting the extreme points in the two outermost bins. The above procedure was repeated by rotating the cluster in increments of 10 degrees around the y-axis to a total of 180 degrees. In the body Gaussian curvature KGK1K2, K1 and K2 denoted the maximum and minimum curvatures. The end-to-end distance Q of the structure was defined as the maximum of the computed Qp values.

2.2.2. Information Entropy of Aggregates (H)

To quantify the changes in the configurational entropy associated with the packing of copolymers within various intermediate assemblies leading to vesiculation, we calculated the probability distribution function (pdf) of the extension of the hydrophobic (PB) and hydrophilic (PEO) segments of the polymer chains. First, the end-to-end distance of each chain (q) in the cluster was calculated. Then, the associated structure parameter for each chain was calculated as sqmaxq, where qmax was the maximum value in q. The pdf of s was normalized such that
1 q m a x q m i n 0 s m p s d s = 1 ,
where sm = qmaxqmin. Given p(s), we calculated the information entropy H of the structure using its standard definition as H(s) ≡ ∑ p(s) ln (p(s)), where the discrete summation was carried out for the ensemble of copolymers within a cluster.

2.2.3. Pair Correlation Function (g(r))

Pair correlation function (g(r)) describes how the density of B, EO, or W beads varies as a function of distance from the reference point of a cluster along the reference direction. A larger value of g signifies tighter packing of the beads. In general, g shows (local density)/(overall density). For spherical structures, the reference point is the center of mass and the reference direction is radial direction. Then, g(r) is given by
g ( r ) = < ρ ( r ) > ρ = 1 ρ i = 1 N δ 1 ( r i r ) δ 2 ( r + r r i ) 1 V ( r )
where < ρ ( r ) > is the average density of beads at a distance r from the center of mass, ρ = 3 N m / ( 4 π r m a x 3 ) is the average cluster density, Nm is the total number of monomer beads contained in the cluster, r m a x is the distance of the farthest bead from the center of mass of the cluster, V ( r ) = ( 4 / 3 ) π ( r + r ) 3 r 3 and
δ 1 ( x ) = 0 ,   if   x < 0 1 ,   if   x 0   and   δ 2 ( x ) = 0 ,   if   x 0 1 ,   if   x > 0 .
Equation (13) means beads at the inner boundary are included while those at the outer boundary are excluded. In the calculation, we used ∆r = 0.1 nm.
For lamellae and cylindrical structures, similar calculation methods were used. The difference was that lamellae and cylindrical structures were divided into small cuboids or disks along the length direction. For each cuboid, the reference point was the lower surface, and the reference direction was thickness direction. For each disk, the reference point was its center, and the reference direction was the radial direction. The g(r) of an entire structure was calculated as the average of the g(r) data obtained for the individual slices.

3. Results and Discussion

3.1. Morphology Evolution

We study the evolution of nearly spherical polymersomes made up of symmetric BCP chains, which consist of six PB and six PEO beads corresponding to a PEO mole fraction of 50%. The BCP concentration is 14.4 wt%. The chain lengths used in the simulations translate to a polymer molecular weight of 590 Da. Hence, these are relatively short polymers with molecular weights appreciably greater than those of surfactants with long hydrocarbon chains, e.g., cetyl-trimethyl ammonium chloride, a C-16 cationic surfactant studied extensively in experiments and our previous CGMD simulations, has a molecular weight of 320 Da.
To quantify the geometry of the structure, we tracked the number of monomers (N), body Gaussian curvature (KG), and end-to-end distance of the self-assembled structure (Q) as functions of time. Representative results are shown in Figure 2 and Figure 3. Vesiculation is facilitated through the following pathway: the formation of well-defined spherical aggregates from a homogeneous copolymer solution, merger of spherical aggregates to form rodlike micelles (t < 100 ns, Figure 2A–F), fusion of rods to form longer wormlike micelles (100 ns < t < 120 ns, Figure 2F–H), flattening of flexible cylindrical structures to form rectangular lamellae (bilayers) that further reorganize into disk-shaped lamellae (120 ns < t < 250 ns, Figure 2H–N), bending and curving of disklike lamellae leading to cavity formation, and closure of the cavity to form vesicles (250 ns < t < 375 ns, Figure 2N–Q). A figure (Figure S6) of this process is provided in Supplementary Materials.
The early stage of vesiculation is characterized by copolymer aggregation into spherical micelles. Each one of the three step changes seen in N in Figure 3 corresponds to a merger event. First, small spherical aggregates consisting of ~1000 monomers coalesce to form a rodlike micelle. Second, two rodlike micelles approach each other at an obtuse angle (see Figure 2F) and combine into a longer wormlike micelle (WM) by an endcap opening mechanism. Specifically, it is energetically favorable for cylindrical micelles to merge as this process results in the loss of the curvature energy associated with their endcaps [72]. Third, the wormlike micelle accommodates a small spherical aggregate and grows larger. The first two mergers occur at terminal positions and are accompanied by an increase in Q while the third merger occurs on the side of the rodlike micelle resulting in little change in Q. For the rodlike structures shown in Figure 2E–H, the positive maximum body curvature appears along the major axis, and the negative minimum body curvature appears along the minor axis. Hence, KG is negative.
After the third merging process, the maximum and standard deviation in N is 0.28% and 0.06%, respectively. Consequently, vesicle formation from the rodlike structure occurs through a series of structure re-organization at a practically constant aggregation number n = N/12, as seen from Figure 3. For 140 ns ≤ t ≤ 180 ns, the reordering of the copolymers within the flexible WM results in the formation of a rectangular lamella (RL). The WM-to-RL transition is accompanied by a continuous decrease in Q. Early in the process, due to the initial angle of approach between the two rodlike micelles during the merging process shown in Figure 2G, a topological twist is formed within the RL, thus creating a saddle point: see row 3 in Table 1. For 180 ns < t ≤ 215 ns, the cluster rearranges to a flat circular disk lamella (DL) to increase the symmetry. During this transition, KG approaches 0, and Q decreases further. For 215 ns < t ≤ 375 ns, spontaneous curvature fluctuations of the DL initiate a curling process. Subsequently. the structure begins to bend and form a spherical cavity (C). This is accompanied by a reduction in Q and an increase in KG. For t > 375 ns, a stable spherical bilayer structure (vesicle, V) is established by the closure of the cavity. Simulations predict that the bending process is characterized by a larger density of monomers (ethylene oxide beads) in the internal leaflet: see Figure 4.

3.2. Prominent Intermediate Structures

Table 1 summarizes the five prominent structures observed along the vesiculation pathway, namely the wormlike micelle (WM), rectangular lamella (RL), disklike lamella (DL), cavity (C), and vesicle (V). Images consisting of blue and red domains depict the entire copolymer, while those with red surfaces reveal only the PB segments. The table also shows the particle pair correlation functions, g(r), and the probability distribution functions, p(s), for the structure parameters corresponding to the segmental stretches of PB (qPB) and PEO (qPEO). Variations in g(r) for PB (hydrophobic), PEO (hydrophilic), and water are shown in red, blue, and black lines, respectively. These functions quantify molecular density profiles within the structures and provide information on the relative sizes of the hydrophilic and hydrophobic domains. A WM is devoid of water inside and has a PB core surrounded by a PEO shell exposed to water. RL, DL, and C all exhibit a PEO-PB-PEO bilayer structure. Vesicles possess a water core and a PEO-PB-PEO bilayer shell. It is intriguing that p(s) variations for all structures are skewed toward lower s (larger segmental stretch q) values and better fitted to a log-normal rather than Gaussian (normal) distributions, discounting the second minor peak seen for PEO. The fitting parameters are provided in Supplementary Materials as Table S4. This second peak in the tail of the pdfs corresponds to ring-like configurations primarily located at the overlapping PEO–PB interfaces (see Figure S1 in Supplementary Materials): see paragraph below. Diffuse (overlapping) interfaces predicted by our simulations are consistent with those observed in previous CGMD simulations of flexible diblock BCPs [4,63]. The average segmental stretch is expectedly greater for the hydrophilic component (PEO). The p(s) of PB features relatively long and low-probability tails, which are a consequence of large chain length fluctuations in the less extended (coiled) PB configurations.
To understand the log normal variations in s, we conducted simulations to track the equilibrium configuration dynamics of a single chain in a solution. An analysis of the data collected for the instantaneous values of qPEO and qPB showed that the chain stretch parameter s follows log-normal distribution for PEO and PB segments as well as the entire diblock chain, as shown in Figure 5. The fitting parameters are provided in Supplementary Materials as Table S5. Consequently, the pdfs of s for assemblies of individual polymer chains are also log-normal. The qualitative agreement of p(s) to a log normal distribution might have implications in terms of entropy maximization. Specifically, a log normal distribution maximizes the information entropy of a random variable with prescribed values of both the mean and variance [80]. The average segmental stretch and configurational fluctuations are bounded by the ring-like and fully stretched states of the segments. The skewness and kurtosis in the distributions appear to depend on the aggregate shape. This suggests that the details of molecular packing within the aggregates are topology dependent. Further investigations are needed to quantify such structure–topology relationships.

3.3. Entropic Motifs

To investigate whether vesicle formation is entropy driven, we calculated the information entropy, H, along the vesiculation pathway within the constant N regime. Figure 6a shows that H values for the PB and PEO segments remain practically unchanged during the vesiculation process, with HPEO and HPB = 4.62 ± 0.10 nat and 5.02 ± 0.02 nat, respectively. Figure 6a also shows that HPB exhibits relatively larger fluctuations. Upon examining qmax and qmin (Figure 6b), we find that, for PEO, qmax and qmin remain practically constant along the vesiculation pathway. However, for PB, qmax shows negligible changes while qmin undergoes relatively large fluctuations. This observation can be attributed to the fact that qmax is related to the fully stretched segment length, which does not change substantially. On the other hand, qmin is associated with the shortest coiled configurations. Based on their chemical structures, the parameter Kθ in the bending potential function for PB is smaller than that for PEO (Kθ = 25 kJ/mol and 50 kJ/mol for PB and PEO, respectively) [75,81]. A lower Kθ value results in greater fluctuations of the coiled configurations. In the absence of systematic statistically significant variations in the information entropy, we conclude that entropic driving forces are not responsible for vesiculation in the system studied.

3.4. Energetic Motifs

We now examine the energetics-based mechanism proposed for lamella-to-vesicle transition [27,34]. Towards this end, we compare the relative change in the surface area ΔS ≡ 1 − Sv/Sd associated with the transition where Sv and Sd denote the surface area exposed to the solvent for the vesicle and disklike lamella, respectively. For an ideal lamella of radius Rd and thickness L, Sd = 2πRd2(1 + L/Rd). If such a lamella bends, curves, and folds to form an ideal vesicle of radius Rv, the resulting surface area Sv = 4π [Rv2 + (RvL)2], where RvRd/2. However, this assumes perfect geometric shapes and that the packing density is conserved during the structure transition process. To assess the impact of such non-idealities on ΔS, we also compute Sd and Sv directly from the simulations. The procedure is illustrated in Figure 7 for the ideal (Figure 7a) and actual (Figure 7b) cases for Rd = 9.9 nm and L = 4.7 nm. Figure 7b,c show that the bilayer thickness L remains unchanged during the lamella-to-vesicle transition, as assumed in classical theoretical models. However, the radius of the resulting vesicle is significantly larger than that predicted by ideal geometric conditions, i.e., simulations yield Rv/Rd = 0.73. This implies that the size of the water core and, consequently, the surface area of the vesicles formed are larger than those predicted by models that use ideal geometric representations of the self-assembled structures. This translates to a smaller value of ΔS ≈ 0.18, predicted by the simulations, vs. ΔS ≈ 0.66, calculated for idealized structures.
Figure 8 shows the change in solvent-accessible surface area (SASA) along the vesiculation pathway. Note that SASA is different from the macroscopic surface area calculated from the linear geometric dimensions shown in Figure 7: specifically, as explained in [30], SASA computations incorporate interfacial area variations on the scale of the monomer bead size. Figure 8a shows that the PEO–water interface accounts for nearly 95% of the cumulative SASA. Moreover, SASA decreases continuously during rod-to-RL and DL-to-vesicle transitions. However, as the RL reorganizes into a DL, SASA is unchanged. This observation is consistent with the changes in the packing density during vesiculation shown in Figure 9. Specifically, during the transition from RL to DL, the packing density and SASA are unaffected. Furthermore, although the SASA associated with PB segments is much smaller compared to that of PEO, the relative change in SASA during rod-to-vesicle transition of the PB segments (≈0.3) is significantly greater than that of PEO segments (≈0.1). This observation supports the edge energy hypothesis [1,13,27,29,31,32,33,34,35] that the unfavorable edge energy of the hydrophobic boundaries is the driving force for vesicle formation. The quantitative differences in SASA between theoretical and simulation predictions can be attributed to the shape imperfections of the aggregates and changes in the packing density. Figure 9 clearly shows that the polymers within the vesicle are more compactly packed as compared to those in the lamella. The reduction in SASA manifests as lower energies of a polymer–solvent interaction, as shown in Figure 10. Specifically, the (non-bonded) interaction energy of either PB or PEO segments with water is seen to become less negative during vesiculation.
Finally, the robustness of the above conclusions was established by performing simulations with different initial configurations of the system: see Figures S2 and S3. Further, the isentropic pathway of vesiculation is qualitatively unchanged with respect to ±10% variations in the bond stretching and bending force parameters: see Figures S4 and S5.

4. Conclusions

CGMD simulations of the real-time evolution of vesicles, starting with an initial condition of randomly distributed diblock copolymer chains in an explicit solvent, have been performed to identify intermediate self-assembled states along the vesiculation pathway. The vesiculation process is facilitated by the rapid formation of spherical aggregates, merger of spherical aggregates to form rodlike assemblies, fusion of rods to form wormlike micelles, flattening of a WM to create a rectangular lamella, reorganization of the RL into a disk-shaped lamella, bending and curving of the DL leading to cavity formation, and closure of the cavity to form a vesicle. An analysis of information entropy shows that vesiculation from a wormlike micelle is an isentropic process. SASA decreases during vesiculation, with a more significant relative surface area reduction occurring in the PB–water interface. These predictions support the edge energy hypothesis that vesiculation is accompanied by a loss in the energy of hydrophobic interactions between solvent and copolymers along the lateral interface of the lamella at the expense of gaining curvature energy. However, simulations show that the rod–lamella–vesicle transition is accompanied by an increase in copolymer packing density. Hence, the change in the surface area accompanying vesiculation predicted by the simulations is significantly lower than theoretical estimates based on idealized geometric representations of lamellar and vesicular morphologies.

Supplementary Materials

The following supporting information can be downloaded at: https://www.mdpi.com/2504-5377/8/1/12/s1, Figure S1. The location of ringlike PEO segments in vesicle. Figure S2. Information entropy vs. time obtained from simulations that start with three different initial copolymer distributions. Figure S3. Number of monomers (N) and end-to-end distance (Q) vs. time obtained from simulations that start with three different initial copolymer distributions. Figure S4. Information entropy vs. time obtained from simulations performed with (a) a 10% reduction and (b) a 10% increase in Kb and Kθ as compared to their values used in the simulations reported in Figure 6. Figure S5. Number of monomers (N) and end-to-end distance (Q) vs. time obtained from simulations performed with (a) a 10% reduction and (b) a 10% increase in Kb and Kθ as compared to their values used in the simulations reported in Figure 3. Figure S6. The vesiculation process. Figure S7. 3D visualization of wormlike micelle intermediate. Figure S8. 3D visualization of rectangular lamella intermediate. Figure S9. 3D visualization of disk lamella intermediate. Figure S10. 3D visualization of cavity micelle intermediate. Figure S11. 3D visualization of vesicle. Table S1. Parameter values for non-bonded interactions. Table S2. Parameter values for bond stretching interactions. Table S3. Parameter values for bond bending interactions. Table S4. The lognormal fitting parameters for the pdf figures in Table 1. Table S5. The lognormal fitting parameters for Figure 5. Section S1. Simulation Details. [75,76,82,83,84] are cited in this section.

Author Contributions

Conceptualization, supervision, project administration, funding acquisition, R.S.; methodology, investigation, resources, writing—original draft preparation, writing—review and editing, S.L. and R.S.; software, validation and formal analysis, data curation, visualization, S.L. All authors have read and agreed to the published version of the manuscript.

Funding

This work was supported by the donors of ACS Petroleum Research Fund under grant #66629-ND9.

Data Availability Statement

Supporting data, figures and animations are provided in Supplementary Materials. Requests for additional data may be sent via email to the authors.

Acknowledgments

This study was supported in part through computational resources provided by Syracuse University. The authors gratefully acknowledge the funding from NSF award ACI-1341006 to support research computing in Syracuse. The visualization figures were made with VMD, which is developed with NIH support by the Theoretical and Computational Biophysics group at the Beckman Institute, University of Illinois at Urbana-Champaign. The authors gratefully acknowledge the developer of Gromacs 2018.1 software/manual.

Conflicts of Interest

The authors declare no conflicts of interest.

Correction Statement

Due to an error in article production, incorrect references were previously listed in the main text. This information has been updated and this change does not affect the scientific content of the article.

References

  1. Bleul, R.; Thiermann, R.; Maskos, M. Techniques to Control Polymersome Size. Macromolecules 2015, 48, 7396–7409. [Google Scholar] [CrossRef]
  2. Mohammadi, M.; Ramezani, M.; Abnous, K.; Alibolandi, M. Biocompatible Polymersomes-Based Cancer Theranostics: Towards Multifunctional Nanomedicine. Int. J. Pharm. 2017, 519, 287–303. [Google Scholar] [CrossRef]
  3. Zhang, X.-Y.; Zhang, P.-Y. Polymersomes in Nanomedicine—A Review. Curr. Nanosci. 2017, 13, 124–129. [Google Scholar] [CrossRef]
  4. Discher, D.E.; Ortiz, V.; Srinivas, G.; Klein, M.L.; Kim, Y.; Christian, D.; Cai, S.; Photos, P.; Ahmed, F. Emerging Applications of Polymersomes in Delivery: From Molecular Dynamics to Shrinkage of Tumors. Prog. Polym. Sci. 2007, 32, 838–857. [Google Scholar] [CrossRef] [PubMed]
  5. Burkett, S.L.; Davis, M.E. Mechanism of Structure Direction In the Synthesis of Pure-Silica Zeolites. 2. Hydrophobic Hydration and Structural Specificity. Chem. Mater. 1995, 7, 1453–1463. [Google Scholar] [CrossRef]
  6. Deng, Y.; Wei, J.; Sun, Z.; Zhao, D. Large-Pore Ordered Mesoporous Materials Templated from Non-Pluronic Amphiphilic Block Copolymers. Chem. Soc. Rev. 2013, 42, 4054–4070. [Google Scholar] [CrossRef] [PubMed]
  7. Van Hest, J.C.; Delnoye, D.A.; Baars, M.W.; van Genderen, M.H.; Meijer, E.W. Polystyrene-Dendrimer Amphiphilic Block Copolymers with a Generation-Dependent Aggregation. Science 1995, 268, 1592–1595. [Google Scholar] [CrossRef] [PubMed]
  8. Zhang, L.; Eisenberg, A. Multiple Morphologies of “Crew-Cut” Aggregates of Polystyrene-b-Poly(Acrylic Acid) Block Copolymers. Science 1995, 268, 1728–1731. [Google Scholar] [CrossRef]
  9. Discher, B.M.; Won, Y.-Y.; Ege, D.S.; Lee, J.C.-M.; Bates, F.S.; Discher, D.E.; Hammer, D.A. Polymersomes: Tough Vesicles Made From Diblock Copolymers. Science 1999, 284, 1143–1146. [Google Scholar] [CrossRef]
  10. Mai, Y.; Eisenberg, A. Self-assembly of Block Copolymers. Chem. Soc. Rev. 2012, 41, 5969. [Google Scholar] [CrossRef]
  11. Jain, S.; Bates, F.S. On the Origins of Morphological Complexity in Block Copolymer Surfactants. Science 2003, 300, 460–464. [Google Scholar] [CrossRef] [PubMed]
  12. Won, Y.-Y.; Brannan, A.K.; Davis, H.T.; Bates, F.S. Cryogenic Transmission Electron Microscopy (Cryo-Tem) of Micelles and Vesicles Formed in Water by Poly(Ethylene Oxide)-Based Block Copolymers. J. Phys. Chem. B 2002, 106, 3354–3364. [Google Scholar] [CrossRef]
  13. Israelachvili, J.N. Intermolecular and Surface Forces; Academic Press: Cambridge, MA, USA, 2011. [Google Scholar]
  14. Karayianni, M.; Pispas, S. Block Copolymer Solution Self-Assembly: Recent Advances, Emerging Trends, and Applications. J. Polym. Sci. 2021, 59, 1874–1898. [Google Scholar] [CrossRef]
  15. Li, S.; Byrne, B.; Welsh, J.E.; Palmer, A.F. Self-Assembled Poly(Butadiene)-b-Poly(Ethylene Oxide) Polymersomes as Paclitaxel Carriers. Biotechnol. Prog. 2007, 23, 278–285. [Google Scholar] [CrossRef] [PubMed]
  16. Allen, C.; Eisenberg, A.; Mrsic, J.; Maysinger, D. PCL-b-PEO Micelles as a Delivery Vehicle for Fk506: Assessment of a Functional Recovery of Crushed Peripheral Nerve. Drug Deliv. 2000, 7, 139–145. [Google Scholar] [CrossRef]
  17. Allen, C.; Yu, Y.; Maysinger, D.; Eisenberg, A. Polycaprolactone-b-Poly(Ethylene Oxide) Block Copolymer Micelles as a Novel Drug Delivery Vehicle for Neurotrophic Agents Fk506 and l-685,818. Bioconjugate Chem. 1998, 9, 564–572. [Google Scholar] [CrossRef]
  18. Boucher-Jacobs, C.; Rabnawaz, M.; Katz, J.S.; Even, R.; Guironnet, D. Encapsulation of catalyst in block copolymer micelles for the polymerization of ethylene in aqueous medium. Nat. Commun. 2018, 9, 841. [Google Scholar] [CrossRef]
  19. Cuomo, F.; Ceglie, A.; De Leonardis, A.; Lopez, F. Polymer Capsules for Enzymatic Catalysis in Confined Environments. Catalysts 2018, 9, 1. [Google Scholar] [CrossRef]
  20. Peters, R.J.; Marguet, M.; Marais, S.; Fraaije, M.W.; van Hest, J.C.; Lecommandoux, S. Cascade Reactions in Multicompartmentalized Polymersomes. Angew. Chem. Int. Ed. 2013, 53, 146–150. [Google Scholar] [CrossRef] [PubMed]
  21. Wilson, D.A.; Nolte, R.J.; van Hest, J.C. Autonomous Movement of Platinum-Loaded Stomatocytes. Nat. Chem. 2012, 4, 268–274. [Google Scholar] [CrossRef] [PubMed]
  22. Marguet, M.; Bonduelle, C.; Lecommandoux, S. Multicompartmentalized Polymeric Systems: Towards Biomimetic Cellular Structure and Function. Chem. Soc. Rev. 2013, 42, 512–529. [Google Scholar] [CrossRef]
  23. Che, H.; van Hest, J.C. Stimuli-Responsive Polymersomes and Nanoreactors. J. Mater. Chem. B 2016, 4, 4632–4647. [Google Scholar] [CrossRef]
  24. Jacobs, M.L.; Boyd, M.A.; Kamat, N.P. Diblock Copolymers Enhance Folding of a Mechanosensitive Membrane Protein during Cell-free Expression. Proc. Natl. Acad. Sci. USA 2019, 116, 4031–4036. [Google Scholar] [CrossRef]
  25. Li, X.; Cooksey, T.J.; Kidd, B.E.; Robertson, M.L.; Madsen, L.A. Mapping Coexistence Phase Diagrams of Block Copolymer Micelles and Free Unimer Chains. Macromolecules 2018, 51, 8127–8135. [Google Scholar] [CrossRef]
  26. Holder, S.W.; Grant, S.C.; Mohammadigoushki, H. Nuclear Magnetic Resonance Diffusometry of Linear and Branched Wormlike Micelles. Langmuir 2021, 37, 3585–3596. [Google Scholar] [CrossRef]
  27. Israelachvili, J.N.; Mitchell, D.J.; Ninham, B.W. Theory of Self-Assembly of Lipid Bilayers and Vesicles. Biochim. Biophys. Acta (BBA)-Biomembr. 1977, 470, 185–201. [Google Scholar] [CrossRef]
  28. Chen, L.; Shen, H.; Eisenberg, A. Kinetics and Mechanism of the Rod-to-Vesicle Transition of Block Copolymer Aggregates in Dilute Solution. J. Phys. Chem. B 1999, 103, 9488–9497. [Google Scholar] [CrossRef]
  29. Antonietti, M.; Förster, S. Vesicles and liposomes: A Self-Assembly Principle Beyond Lipids. Adv. Mater. 2003, 15, 1323–1333. [Google Scholar] [CrossRef]
  30. Liu, S.; Sureshkumar, R. Morphological Diversity in Diblock Copolymer Solutions: A Molecular Dynamics Study. Colloids Interfaces 2023, 7, 40. [Google Scholar] [CrossRef]
  31. Huang, C.; Quinn, D.; Sadovsky, Y.; Suresh, S.; Hsia, K.J. Formation and Size Distribution of Self-Assembled Vesicles. Proc. Natl. Acad. Sci. USA 2017, 114, 2910–2915. [Google Scholar] [CrossRef] [PubMed]
  32. Srinivas, G.; Discher, D.E.; Klein, M.L. Self-Assembly and Properties of Diblock Copolymers by Coarse-Grain Molecular Dynamics. Nat. Mater. 2004, 3, 638–644. [Google Scholar] [CrossRef] [PubMed]
  33. Srinivas, G.; Shelley, J.C.; Nielsen, S.O.; Discher, D.E.; Klein, M.L. Simulation of Diblock Copolymer Self-Assembly, Using a Coarse-Grain Model. J. Phys. Chem. B 2004, 108, 8153–8160. [Google Scholar] [CrossRef]
  34. Lipowsky, R. The Conformation of Membranes. Nature 1991, 349, 475–481. [Google Scholar] [CrossRef]
  35. Fromherz, P.; Röcker, C.; Rüppel, D. From discoid micelles to spherical vesicles. The concept of edge activity. Faraday Discuss. Chem. Soc. 1986, 81, 39–48. [Google Scholar] [CrossRef]
  36. Noguchi, H.; Gompper, G. Dynamics of vesicle self-assembly and dissolution. J. Chem. Phys. 2006, 125, 164908. [Google Scholar] [CrossRef]
  37. Yuan, H.; Huang, C.; Li, J.; Lykotrafitis, G.; Zhang, S. One-particle-thick, solvent-free, coarse-grained model for biological and biomimetic fluid membranes. Phys. Rev. E 2010, 82, 011905. [Google Scholar] [CrossRef]
  38. Canham, P.B. The Minimum Energy of Bending as a Possible Explanation of the Biconcave Shape of the Human Red Blood Cell. J. Theor. Biol. 1970, 26, 61–81. [Google Scholar] [CrossRef]
  39. Helfrich, W. Elastic Properties of Lipid Bilayers: Theory and Possible Experiments. Z. Naturforschung C 1973, 28, 693–703. [Google Scholar] [CrossRef] [PubMed]
  40. Faucon, J.F.; Mitov, M.D.; Méléard, P.; Bivas, I.; Bothorel, P. Bending Elasticity and Thermal Fluctuations of Lipid Membranes-Theoretical and Experimental Requirements. J. Phys. 1989, 50, 2389–2414. [Google Scholar] [CrossRef]
  41. Evans, E.; Needham, D. Physical Properties of Surfactant Bilayer Membranes: Thermal Transitions, Elasticity, Rigidity, Cohesion and Colloidal Interactions. J. Phys. Chem. 1987, 91, 4219–4228. [Google Scholar] [CrossRef]
  42. Cooke, I.R.; Deserno, M. Coupling Between Lipid Shape and Membrane Curvature. Biophys. J. 2006, 91, 487–495. [Google Scholar] [CrossRef]
  43. Fogel, A.L.; Ravichandran, A.; Mani, S.; Upadhyay, B.; Khare, R.; Morgan, S.E. Water structure and mobility in acrylamide copolymer glycohydrogels with galactose and Siloxane Pendant Groups. J. Polym. Sci. Part B Polym. Phys. 2019, 57, 584–597. [Google Scholar] [CrossRef]
  44. Walter, A.; Vinson, P.K.; Kaplun, A.; Talmon, Y. Intermediate structures in the cholate-phosphatidylcholine vesicle-micelle transition. Biophys. J. 1991, 60, 1315–1325. [Google Scholar] [CrossRef]
  45. Vinson, P.K.; Talmon, Y.; Walter, A. Vesicle-micelle transition of phosphatidylcholine and octyl glucoside elucidated by cryo-transmission electron microscopy. Biophys. J. 1989, 56, 669–681. [Google Scholar] [CrossRef]
  46. Davies, T.S.; Ketner, A.M.; Raghavan, S.R. Self-assembly of surfactant vesicles that transform into viscoelastic wormlike micelles upon heating. J. Am. Chem. Soc. 2006, 128, 6669–6675. [Google Scholar] [CrossRef]
  47. Markvoort, A.J.; van Santen, R.A.; Hilbers, P.A. Vesicle shapes from molecular dynamics simulations. J. Phys. Chem. B 2006, 110, 22780–22785. [Google Scholar] [CrossRef]
  48. Markvoort, A.J.; Pieterse, K.; Steijaert, M.N.; Spijker, P.; Hilbers, P.A. The Bilayer−Vesicle Transition is Entropy Driven. J. Phys. Chem. B 2005, 109, 22649–22654. [Google Scholar] [CrossRef] [PubMed]
  49. He, X.; Schmid, F. Spontaneous formation of complex micelles from a homogeneous solution. Phys. Rev. Lett. 2008, 100, 137802. [Google Scholar] [CrossRef] [PubMed]
  50. He, X.; Schmid, F. Dynamics of spontaneous vesicle formation in dilute solutions of Amphiphilic Diblock Copolymers. Macromolecules 2006, 39, 2654–2662. [Google Scholar] [CrossRef]
  51. Sevink, G.J.; Zvelindovsky, A.V. Self-assembly of complex vesicles. Macromolecules 2005, 38, 7502–7513. [Google Scholar] [CrossRef]
  52. Ye, X.; Khomami, B. Self-Assembly of Linear Diblock Copolymers in Selective Solvents: From Single Micelles to Particles with Tri-Continuous Inner Structures. Soft Matter 2020, 16, 6056–6062. [Google Scholar] [CrossRef] [PubMed]
  53. Li, Z.; Dormidontova, E.E. Equilibrium Chain Exchange Kinetics in Block Copolymer Micelle Solutions by Dissipative Particle Dynamics Simulations. Soft Matter 2011, 7, 4179. [Google Scholar] [CrossRef]
  54. Javan Nikkhah, S.; Turunen, E.; Lepo, A.; Ala-Nissila, T.; Sammalkorpi, M. Multicore assemblies from three-component linear homo-copolymer systems: A coarse-grained modeling study. Polymers 2021, 13, 2193. [Google Scholar] [CrossRef] [PubMed]
  55. Ortiz, V.; Nielsen, S.O.; Discher, D.E.; Klein, M.L.; Lipowsky, R.; Shillcock, J. Dissipative particle dynamics simulations of polymersomes. J. Phys. Chem. B 2005, 109, 17708–17714. [Google Scholar] [CrossRef] [PubMed]
  56. Shillcock, J.C. Spontaneous vesicle self-assembly: A mesoscopic view of membrane dynamics. Langmuir 2012, 28, 541–547. [Google Scholar] [CrossRef] [PubMed]
  57. Xiao, M.; Liu, J.; Yang, J.; Wang, R.; Xie, D. Biomimetic membrane control of block copolymer vesicles with tunable wall thickness. Soft Matter 2013, 9, 2434. [Google Scholar] [CrossRef]
  58. Li, Y.; Zhang, H.; Wang, Z.; Bao, M. Micelle-vesicle transitions in catanionic mixtures of SDS/DTAB induced by salt, temperature, and selective solvents: A dissipative particle dynamics simulation study. Colloid Polym. Sci. 2014, 292, 2349–2360. [Google Scholar] [CrossRef]
  59. Luo, Z.; Li, Y.; Wang, B.; Jiang, J. Ph-sensitive vesicles formed by amphiphilic grafted copolymers with tunable membrane permeability for drug loading/release: A multiscale simulation study. Macromolecules 2016, 49, 6084–6094. [Google Scholar] [CrossRef]
  60. Wu, S.; Lu, T.; Guo, H. Dissipative particle dynamic simulation study of lipid membrane. Front. Chem. China 2010, 5, 288–298. [Google Scholar] [CrossRef]
  61. Feng, X.; Yan, N.; Jin, J.; Jiang, W. Disassembly of amphiphilic AB block copolymer vesicles in selective solvents: A molecular dynamics simulation study. Macromolecules 2023, 56, 2560–2567. [Google Scholar] [CrossRef]
  62. Lyubimov, I.; Wessels, M.G.; Jayaraman, A. Molecular dynamics simulation and prism theory study of assembly in solutions of amphiphilic bottlebrush block copolymers. Macromolecules 2018, 51, 7586–7599. [Google Scholar] [CrossRef]
  63. Chakraborty, K.; Shinoda, W.; Loverde, S.M. Molecular Simulation of the Shape Deformation of a Polymersome. Soft Matter 2020, 16, 3234–3244. [Google Scholar] [CrossRef]
  64. Sun, X.; Pei, S.; Wang, J.; Wang, P.; Liu, Z.; Zhang, J. Coarse-grained molecular dynamics simulation study on spherical and tube-like vesicles formed by amphiphilic copolymers. J. Polym. Sci. Part B Polym. Phys. 2017, 55, 1220–1226. [Google Scholar] [CrossRef]
  65. Marrink, S.J.; Mark, A.E. Molecular dynamics simulation of the formation, structure, and dynamics of small phospholipid vesicles. J. Am. Chem. Soc. 2003, 125, 15233–15242. [Google Scholar] [CrossRef]
  66. Wu, R.; Deng, M.; Kong, B.; Yang, X. Coarse-grained molecular dynamics simulation of ammonium surfactant self-assemblies: Micelles and vesicles. J. Phys. Chem. B 2009, 113, 15010–15016. [Google Scholar] [CrossRef]
  67. Sambasivam, A.; Dhakal, S.; Sureshkumar, R. Structure and Rheology of Self-Assembled Aqueous Suspensions of Nanoparticles and Wormlike Micelles. Mol. Simul. 2017, 44, 485–493. [Google Scholar] [CrossRef]
  68. Sambasivam, A.; Sangwai, A.V.; Sureshkumar, R. Self-Assembly of Nanoparticle–Surfactant Complexes with Rodlike Micelles: A Molecular Dynamics Study. Langmuir 2016, 32, 1214–1219. [Google Scholar] [CrossRef]
  69. Sambasivam, A.; Sangwai, A.V.; Sureshkumar, R. Dynamics and Scission of Rodlike Cationic Surfactant Micelles in Shear Flow. Phys. Rev. Lett. 2015, 114, 8302. [Google Scholar] [CrossRef]
  70. Dhakal, S.; Sureshkumar, R. Anomalous Diffusion and Stress Relaxation in Surfactant Micelles. Phys. Rev. E 2017, 96, 2605. [Google Scholar] [CrossRef]
  71. Dhakal, S.; Sureshkumar, R. Uniaxial Extension of Surfactant Micelles: Counterion Mediated Chain Stiffening and a Mechanism of Rupture by Flow-Induced Energy Redistribution. ACS Macro Lett. 2015, 5, 108–111. [Google Scholar] [CrossRef]
  72. Dhakal, S.; Sureshkumar, R. Topology, Length Scales, and Energetics of Surfactant Micelles. J. Chem. Phys. 2015, 143, 024905. [Google Scholar] [CrossRef] [PubMed]
  73. Sangwai, A.V.; Sureshkumar, R. Binary Interactions and Salt-Induced Coalescence of Spherical Micelles of Cationic Surfactants From Molecular Dynamics Simulations. Langmuir 2011, 28, 1127–1135. [Google Scholar] [CrossRef] [PubMed]
  74. Sangwai, A.V.; Sureshkumar, R. Coarse-Grained Molecular Dynamics Simulations of the Sphere to Rod Transition in Surfactant Micelles. Langmuir 2011, 27, 6628–6638. [Google Scholar] [CrossRef] [PubMed]
  75. Marrink, S.J.; Risselada, H.J.; Yefimov, S.; Tieleman, D.P.; de Vries, A.H. The Martini Force Field:  Coarse Grained Model for Biomolecular Simulations. J. Phys. Chem. B 2007, 111, 7812–7824. [Google Scholar] [CrossRef]
  76. Darré, L.; Machado, M.R.; Pantano, S. Coarse-grained models of water. WIREs Comput. Mol. Sci. 2012, 2, 921–930. [Google Scholar] [CrossRef]
  77. Martínez, L.; Andrade, R.; Birgin, E.G.; Martínez, J.M. Packmol: A package for building initial configurations for molecular dynamics simulations. J. Comput. Chem. 2009, 30, 2157–2164. [Google Scholar] [CrossRef]
  78. Eisenhaber, F.; Lijnzaad, P.; Argos, P.; Sander, C.; Scharf, M. The double cubic lattice method: Efficient approaches to numerical integration of surface area and volume and to Dot surface contouring of Molecular Assemblies. J. Comput. Chem. 1995, 16, 273–284. [Google Scholar] [CrossRef]
  79. Pratt, V. Direct least-squares fitting of algebraic surfaces. ACM SIGGRAPH Comput. Graph. 1987, 21, 145–152. [Google Scholar] [CrossRef]
  80. Park, S.Y.; Bera, A.K. Maximum Entropy Autoregressive conditional heteroskedasticity model. J. Econom. 2009, 150, 219–230. [Google Scholar] [CrossRef]
  81. Yang, Y. Structure, Dynamics and Rheology of Polymer Solutions from Coarse-Grained Molecular Dynamics: Effects of Polymer Concentration, Solvent Quality and Geometric Confinement. Ph.D. Dissertation, Syracuse University, Syracuse, NY, USA, 2015. [Google Scholar]
  82. Lee, H.; de Vries, A.H.; Marrink, S.-J.; Pastor, R.W. A coarse-grained model for polyethylene oxide and polyethylene glycol: Conformation and hydrodynamics. J. Phys. Chem. B 2009, 113, 13186–13194. [Google Scholar] [CrossRef]
  83. Bussi, G.; Donadio, D.; Parrinello, M. Canonical Sampling Through Velocity Rescaling. J. Chem. Phys. 2007, 126, 014101. [Google Scholar] [CrossRef] [PubMed]
  84. Berendsen, H.J.; Postma, J.P.; van Gunsteren, W.F.; DiNola, A.; Haak, J.R. Molecular Dynamics with Coupling to An External Bath. J. Chem. Phys. 1984, 81, 3684–3690. [Google Scholar] [CrossRef]
Figure 1. Illustration of the algorithm used to compute the body curvature KG and end-to-end distance Q for a disklike lamella (1(ac)) and a cavity (1(df)). (a,d) show 3−d mesh representations of the structures constructed using the position vectors of the outermost monomers. (b,e) show the cross sections marked in (a,d), respectively. The blue and red dots represent the monomer beads and the center of mass of the beads within the intervals marked by the solid black lines in (b,e), respectively. (e,f) show the best fit circle with radius equal to the inverse curvature for the lamella and cavity, respectively.
Figure 1. Illustration of the algorithm used to compute the body curvature KG and end-to-end distance Q for a disklike lamella (1(ac)) and a cavity (1(df)). (a,d) show 3−d mesh representations of the structures constructed using the position vectors of the outermost monomers. (b,e) show the cross sections marked in (a,d), respectively. The blue and red dots represent the monomer beads and the center of mass of the beads within the intervals marked by the solid black lines in (b,e), respectively. (e,f) show the best fit circle with radius equal to the inverse curvature for the lamella and cavity, respectively.
Colloids 08 00012 g001
Figure 2. Intermediate states during vesiculation. The blue and red beads represent PEO and PB monomers, respectively. (A): homogeneous solution; (BE): formation of spherical and cylindrical micelles; (EH): merger of spherical and cylindrical structures to form a flexible wormlike micelle (WM); (HL): reorganization of WM to form a rectangular lamella (RL); (LN): formation of disklike lamella (DL); (OQ): shape fluctuations of DL leading cavity and vesicle formation. For states (O,P), top and side views are shown.
Figure 2. Intermediate states during vesiculation. The blue and red beads represent PEO and PB monomers, respectively. (A): homogeneous solution; (BE): formation of spherical and cylindrical micelles; (EH): merger of spherical and cylindrical structures to form a flexible wormlike micelle (WM); (HL): reorganization of WM to form a rectangular lamella (RL); (LN): formation of disklike lamella (DL); (OQ): shape fluctuations of DL leading cavity and vesicle formation. For states (O,P), top and side views are shown.
Colloids 08 00012 g002
Figure 3. Number of monomers (N), body Gaussian curvature (KG), and end-to-end distance (Q) vs. time. Blue and red colors represent PEO and PB domains. For the discoid lamella, top and side views are shown. KG = 0 along the dotted horizontal line. The first two data points of Q (the stage that has many micelles) are for the largest precursor micelle that forms a well-defined vesicle.
Figure 3. Number of monomers (N), body Gaussian curvature (KG), and end-to-end distance (Q) vs. time. Blue and red colors represent PEO and PB domains. For the discoid lamella, top and side views are shown. KG = 0 along the dotted horizontal line. The first two data points of Q (the stage that has many micelles) are for the largest precursor micelle that forms a well-defined vesicle.
Colloids 08 00012 g003
Figure 4. Density of ethylene oxide beads in the top and bottom leaflets of the bilayer as a function of time. The dotted lines are moving averages.
Figure 4. Density of ethylene oxide beads in the top and bottom leaflets of the bilayer as a function of time. The dotted lines are moving averages.
Colloids 08 00012 g004
Figure 5. Probability distribution function of the segmental stretch parameter for equilibrium (a). PB segment, (b). PEO segment, and (c). PB-PEO chain. Typical configurations at different locations in the pdf are shown.
Figure 5. Probability distribution function of the segmental stretch parameter for equilibrium (a). PB segment, (b). PEO segment, and (c). PB-PEO chain. Typical configurations at different locations in the pdf are shown.
Colloids 08 00012 g005
Figure 6. (a) Information entropy and (b) maximum and minimum segmental extension vs. time.
Figure 6. (a) Information entropy and (b) maximum and minimum segmental extension vs. time.
Colloids 08 00012 g006
Figure 7. (a) Schematic diagram of disk-to-vesicle transition. (b) Disklike lamella and (c) vesicle from simulation.
Figure 7. (a) Schematic diagram of disk-to-vesicle transition. (b) Disklike lamella and (c) vesicle from simulation.
Colloids 08 00012 g007
Figure 8. (a) SASA and (b) relative SASA vs. time. The red line is guide to the eye.
Figure 8. (a) SASA and (b) relative SASA vs. time. The red line is guide to the eye.
Colloids 08 00012 g008
Figure 9. Packing density vs. time. The red line is guide to the eye.
Figure 9. Packing density vs. time. The red line is guide to the eye.
Colloids 08 00012 g009
Figure 10. Total non-bonded pairwise interaction energy for (a) water-PB (b) water-PEO. Black squares denote the simulation data and red dotted line represent the moving average of the data.
Figure 10. Total non-bonded pairwise interaction energy for (a) water-PB (b) water-PEO. Black squares denote the simulation data and red dotted line represent the moving average of the data.
Colloids 08 00012 g010
Table 1. Quantitative descriptors of various morphologies. In the particle pair correlation functions, black, red, and blue lines represent water, PB, and PEO, respectively. See Figures S7–S11 for corresponding 3D visualizations.
Table 1. Quantitative descriptors of various morphologies. In the particle pair correlation functions, black, red, and blue lines represent water, PB, and PEO, respectively. See Figures S7–S11 for corresponding 3D visualizations.
StructurePair Correlation FunctionProbability Distribution Function, PEOProbability Distribution Function, PB
Wormlike micelle
Colloids 08 00012 i016
Colloids 08 00012 i001Colloids 08 00012 i002Colloids 08 00012 i003
Rectangular lamella
Colloids 08 00012 i017
Colloids 08 00012 i004Colloids 08 00012 i005Colloids 08 00012 i006
Disk lamella
Colloids 08 00012 i018
Colloids 08 00012 i007Colloids 08 00012 i008Colloids 08 00012 i009
Cavity micelle
Colloids 08 00012 i019
Colloids 08 00012 i010Colloids 08 00012 i011Colloids 08 00012 i012
Vesicle
Colloids 08 00012 i020
Colloids 08 00012 i013Colloids 08 00012 i014Colloids 08 00012 i015
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

Liu, S.; Sureshkumar, R. Energetic and Entropic Motifs in Vesicle Morphogenesis in Amphiphilic Diblock Copolymer Solutions. Colloids Interfaces 2024, 8, 12. https://doi.org/10.3390/colloids8010012

AMA Style

Liu S, Sureshkumar R. Energetic and Entropic Motifs in Vesicle Morphogenesis in Amphiphilic Diblock Copolymer Solutions. Colloids and Interfaces. 2024; 8(1):12. https://doi.org/10.3390/colloids8010012

Chicago/Turabian Style

Liu, Senyuan, and Radhakrishna Sureshkumar. 2024. "Energetic and Entropic Motifs in Vesicle Morphogenesis in Amphiphilic Diblock Copolymer Solutions" Colloids and Interfaces 8, no. 1: 12. https://doi.org/10.3390/colloids8010012

APA Style

Liu, S., & Sureshkumar, R. (2024). Energetic and Entropic Motifs in Vesicle Morphogenesis in Amphiphilic Diblock Copolymer Solutions. Colloids and Interfaces, 8(1), 12. https://doi.org/10.3390/colloids8010012

Article Metrics

Back to TopTop