MOLECULAR THERMODYNAMICS OF MICELLIZATION : MICELLE SIZE DISTRIBUTIONS AND GEOMETRY TRANSITIONS

Surfactants are amphiphilic molecules that can spontaneously self-assemble in solution, forming structures known as micelles. Variations in temperature, pH, and electrolyte concentration imply changes in the interactions between surfactants and micelle stability conditions, including micelle size distribution and micelle shape. Here, molecular thermodynamics is used to describe and predict conditions of micelle formation in surfactant solutions by directly calculating the minimum Gibbs free energy of the system, corresponding to the most stable condition of the surfactant solution. In order to find it, the proposed methodology takes into account the micelle size distribution and two possible geometries (spherical and spherocylindrical). We propose a numerical optimization methodology where the minimum free energy can be reached faster and in a more reliable way. The proposed models predict the critical micelle concentration well when compared to experimental data, and also predict the effect of salt on micelle geometry transitions.


INTRODUCTION
Self-assembly is a broad term applied to spontaneous organization of substances under appropriate conditions and proportions.It is a reversible process and represents a condition of thermodynamic equilibrium, which is relevant to a large variety of phenomena and systems, such as crystal formation, colloidal systems, lipid bilayers, among others (Whitesides and Boncheva, 2002).When surfactants self-assemble, the structures formed are called micelles.Micellization phenomena have been of great interest to both the academic community and industry.It is involved in several applications and processes.We can mention its application in personal care products, paints, processed food, and separations in the petrochemical industry.Intensification of studies in this field is reinforced by the accessibility of critical micelle concentration data and its great importance for the physicochemical behavior of surfactant solutions.Furthermore, several thermodynamic and transport properties of surfactant solutions are affected by the size and shape of micelles (Iyer and Blankschtein, 2012).For example, the solubility of micelles increase with their size (Rusanov, 2014), and micelle shape transitions significantly affect the viscosity of the solution (Kamranfar and Jamialahmadi, 2014).All these facts contribute to the development of models that can predict the behavior of surfactant solutions and also can be used as a tool for designing new surfactants.

Brazilian Journal of Chemical Engineering
The molecular thermodynamics approach to describe the micellization phenomena was first proposed by Tanford (Tanford, 1974).Tanford analyzed the adequacy of the Gibbs thermodynamic equilibrium to describe the self-assembly of amphiphilic molecules.In Tanford's work, the equilibrium state is calculated by imposing the necessary condition of a stoichiometric linear relation of chemical potential of all species.That implies that the sum of the Gibbs free energies of reactants (amphiphilic molecules) must be the equal to the sum of the Gibbs free energies of products (micelles).
Molecular thermodynamic models aim to be completely predictive.They describe the micellization phenomena only based on information about the molecular characteristics of substances and conditions of the medium, such as, for example, temperature, surfactant concentration, and ionic strength (Goldsipe and Blankschtein, 2007).Nagarajan and collaborators made important contributions to the development of molecular thermodynamic models used to describe the self-assembly of surfactant solutions.In their models (Nagarajan and Ruckenstein, 1991;Nagarajan, 1993;Nagarajan, 2002) they analyzed small and monodisperse micelles and also the transition between two micelle geometries.However, in those works, only the necessary condition for equilibrium is considered.The Nagarajan model was improved by Moreira and Firoozabadi (2009) when they started considering the minimum of Gibbs free energy as a condition for the thermodynamic equilibrium of the solution.However, Moreira and Firoozabadi (2009) used the maximum term approach, which assumes that the micelle size distribution can be represented by one characteristic micelle size and a single number of micelles formed.Therefore, Moreira and Firoozabadi (2009) did not take into account the fact that micelles can be formed with a large distribution of sizes.They also assumed that only small micelles could be formed, with either a spherical or globular shape.In their subsequent works (Moreira and Firoozabadi, 2010;Moreira and Firoozabadi, 2012;Lukanov and Firoozabadi, 2014) those aspects were kept the same.Therefore, to the best of our knowledge, the present work is the only one that considers the minimum Gibbs free energy as the most stable state together with the micelles' size distributions and also analyzes the transition between small (spherical) and rod-like (spherocylindrical) micelles as a function of temperature and salt concentration.
In this work we present an approach where the minimum of total Gibbs free energy is obtained for aqueous surfactant solutions (with and without electrolytes) considering the micelle size distributions and the more stable geometry (spherical and spherocylindrical) of micelles.Here, the stable state for a given set of temperature (T), pressure (P) and composition (N), obtained by the minimization of total Gibbs free energy of the system, gives us information about the most stable micelle geometry, micelles sizes and the critical micelle concentration (CMC) of the solution.

MODEL
The model used to calculate the total Gibbs free energy of the system is similar to the one proposed by Nagarajan and Ruckenstein (1991) and Moreira and Firoozabadi (2009).Given a global specification of T (temperature), P (pressure), SA N (total number of surfactant molecules in the solution) and w N (number of water molecules), the Gibbs free energy is calculated as a sum of two contributions: free energy of formation and free energy of mixture. where  and o g  are the standard chemical potentials of water, free surfactant and micelle with aggregation number g, Ng is the number of micelles, k is the Boltzmann constant, and Xw, X1A and Xg are the mole fractions of water, free surfactant and micelles.It is possible to reorganize the previous expression by separating the terms that depend only on fixed variables: T, P, NSA, and Nw, and dividing the expression by kT: Equation ( 2) represents the expression to be minimized to predict the most stable state of the system.(1991) and Moreira and Firoozabadi (2009).

The free energy of micellization
When the micelles size distribution is narrow enough, it is common to use the Maximum Term approximation.when the Maximum Term approximation is used, Equation (2) becomes:

METHODOLOGY Critical Micelle Concentration (CMC)
Experimentally, the critical micelle concentration (CMC) is defined as the concentration of surfactant where a sharp change in any property of the solution is observed.This information is available in the literature for a wide variety of surfactants.Therefore we decided to validate our model by comparing our predictive calculations with experimental CMC data.To obtain the CMC, we perform different simulations of a solution containing an increasing amount of surfactant.Then, for each solution with a different amount of surfactant added we performed the minimization of the Gibbs free energy.This procedure permits one to relate the total number of surfactant molecules added (NSA) with the number of free surfactant molecules in the solution (N1A), as is graphically presented in Figure 1.The CMC is then defined as the concentration of surfactant added (NSA) where an inflexion of the curve is observed.To obtain this point automatically, a regularization function is used which relates N1A as a function of NSA, according to Equation (4).
where NSA * , determined by non-linear regression, corresponds to the critical number of surfactant molecules added to the solution under the conditions T and P, and ξ is the regularization parameter, where 1 0    .The parameters a, b, and c are obtained from the fitting of the data to the regularization function.
To perform the parameter fitting we used the deterministic method Sequential Quadratic Programming (SQP).After obtaining the concentration of surfactant added (NSA), we can calculate the critical micelle concentration using the number of water molecules assumed in the simulation (NW).

Micelle Size Distribution
To develop the methodology for the micelle size distributions, first we verify what kind of behavior this distribution can have.For that, we consider a set of optimization variables, where it is assumed that we could obtain micelles with an aggregation number from 2 to 100.To construct the next figure, we define as optimization variables the vector of the number of micelles formed for each aggregation number g.Thus, for this specific case we did not pre-set any shape for the micelle size distribution.Figure 2 shows the tendency obtained for the micelle size distribution for an aqueous solution of sodium dodecyl sulfate (SDS) at a concentration ten times higher than its critical micelle concentration.It is interesting to observe that the problem converges to a micelle size distribution with a behavior similar to a Gaussian distribution.
From now on in this work, a Gaussian shape is assumed for all the micelle size distributions.Therefore, it is possible to define the total number of micelles formed, the mean (or expectation) aggregation number of the distribution and its standard deviation as important variables for calculating the minimum free energy of the system.
To validate our methodology for the micelle size distribution, first we compared it with the maximum term approach, both considering the minimum Gibbs free energy to define the most stable state.Because we calculate the critical micelle concentration, only small micelles are considered (spherical and globular shapes only).The geometrical relations used for this section are similar to those presented by Nagarajan and Ruckenstein (1991).Using our methodology, three optimization variables, as mentioned before, are used to calculate each minimum of Gibbs energy.Using the Maximum Term Approximation, we have two optimization variables: the number of micelles formed, Ng, and the specific aggregation number, g.

Geometry Optimization
As in all optimization problems, we must pay attention to the correct definition of the optimization variables to be used in the numerical algorithm.The incorrect selection of variables, if they do not satisfy some constraints or correspond to an over specification, may lead to erroneous results and misinterpretations.It is even more important when using nondeterministic algorithms because the analyticity of the functions and equality constraints are not taken into account.This can result in objective functions with several local minima.
Analyzing the geometric relations presented in Table 1, we can observe that, if the length of the cylindrical part of the micelle (Lc) is zero, the geometry reduces to a spherical micelle.For cases where the geometric transition is considered, we define as optimization variables the number of micelles formed (Ng), the number of surfactant molecules in the cylindrical part of the micelle (gcyl), the radius of the spherical ends of the micelles (Rs), and the length of the cylindrical part (Lc).The number of surfactant molecules in the spherical ends of the micelle is then obtained from the geometric relations presented in Table 1, which assumes the form of: To reduce the optimization search region, we defined limits for the optimization variables.The maximum value of the radius, Rs, is assumed to be equal to the extended length of the surfactant tail, ensuring that there is no empty space inside the micelle.Besides that, Rs must be equal to or greater than the radius of the cylindrical part of the micelle to ensure that the amount of surfactant molecules in the spherical ends is positive.Figure 3 presents a surface of the Gibbs energy, to be minimized, fixing values for Lc and Rs.We can observe that the function has a smooth behavior and a well-defined global minimum.

Numerical Method
As presented in the previous section, calculating the critical micelle concentration involves repeating the Gibbs energy minimization several times for different surfactant concentrations.To reduce the computational time to perform these calculations, we propose the following optimization procedure.First, we define an amount of water to be considered and the temperature of the system.After that, we start our simulations with a small amount of surfactant (to guarantee that we are below the critical micelle concentration).For this first simulation we use the Particle Swarm Optimization method (PSO).The result is then used as an initial guess for the Sequential Quadratic Programming (SQP) deterministic method.Then, we make a small increment in the amount of surfactant and perform the optimization with the deterministic method.The initial guess now, in this case, is the optimum value of the optimization at the lower surfactant concentration.With this methodology, it is only necessary to use the nondeterministic method once, which largely reduces the computational effort.All simulations were performed with the software Matlab® version R2008a.More details about the PSO method can be found elsewhere (Kennedy and Eberhart, 1995) and for the SQP method in Byrd et al. (2000).
A R  Vg is the volume of the micelle, Ag is its surface area, Agδ is its area at a distance δ from the surface, Pf is the packing factor, and R is the radius of the micelle (s for the spherical part and c for the cylindrical part).Variable H is a geometrical parameter.The expression of the volume of the surfactant tail (υs), as well as molecular parameters for different surfactants are presented elsewhere (Nagarajan and Ruckstein, 1991).

RESULTS AND DISCUSSION
Figure 4 presents the critical micelle concentration (CMC) for different surfactants as a function of their tail length.To perform these calculations we assumed that micelles are globular or spherically shaped.The performance of the maximum term method and the micelle size distribution method were compared.We observe that both methodologies are in very good agreement and predict the critical micelle concentrations well.We expect that at low surfactant concentration (close to the CMC) the micelle-size distribution tends to be narrower, compared with high surfactant concentration.For example, in Figure 5 we present the micelle size distribution for different surfactants at different concentrations.We observe relatively wide size distributions for all the cases presented in Figure 5, meaning that, for those cases, the Maximum Term method certainly is not the best approach.We expect wider distributions for high surfactant concentration or for microemulsion systems, in which case the methodology proposed here becomes more adequate to describe those systems.Another aspect that motivates the use of the size distribution approach proposed here is because sizes (mean, variance and other moments of the distribution) are important information for the stability of emulsions.For example, when an electric field is applied to an emulsion, its stability is a function of the diameter of the micelles to the fourth power.
We also proposed a methodology to describe the micelle geometry transitions.The results obtained are presented in Table 2.The surfactant dodecyl trimethylammonium bromide (C12TAB) forms spherical micelles even at high surfactant concentration.However, sodium dodecyl sulfate (SDS) presents a transition from spherical to spherocylindrical as observed experimentally (Victorov and Koroleva, 2014).Another result obtained is that the addition of electrolytes to the surfactant solution causes the formation of spherocylindrical micelles for tetradecyl trimethylammonium bromide (C14TAB), and hexadecyl trimethylammonium bromide (C16TAB).For the nonionic surfactant, decyl(dimethyl)phosphine oxide (C10PO), the model predicts the formation of spherocylindrical micelles even without electrolytes added to the solution, but at higher surfactant concentrations.In general, we are in good qualitative agreement with the trends observed by experiments.

Figure 1 :
Figure 1: Determination of the CMC.Each point represents one minimization process of the Gibbs free energy and the continuous line is the adjusted regulation function.CMC is obtained by parameter fitting of the regularization function represented by Equation (4).

Figure 2 :
Figure 2: Micelle size distribution for a solution of sodium dodecyl sulfate (SDS) at a concentration ten times its critical micelle concentration.Xg is the mole fraction of micelles with the aggregation number g. T = 25 ºC, P = 1 atm, and w N = 3400000.The micelle size distribution converges to a Gaussian-like distribution after the optimization procedure.

Figure 3 :
Figure 3: Typical example of the Gibbs energy, / G kT , as a function of the number of micelles with aggregation number g, Ng, and the number of surfactant molecules in the cylindrical part of the micelle, gcyl.We fixed the values of the length of the cylindrical part of the micelle, Lc, and the radius of the spherical ends of the micelles, Rs, to be equal to 2.89 nm and 1.23 nm, respectively.

Figure 5 :a
Figure 5: Micelle size distributions for different surfactants obtained by the proposed methodology.is the molar fraction of micelles with aggregation number : (a) 16.34 mM SDS solution; (b) 81.7 mM C12TAB solution; (c) 32.68 mM C14TAB solution; (d) 32.68 mM C16TAB solution.We defined = 25 ºC, P = 1 atm, and w N = 3400000.g X g T