Minimum Free-Energy Paths for the Self-Organization of Polymer Brushes

A methodology to calculate minimum free-energy paths based on the combination of a molecular theory and the improved string method is introduced and applied to study the self-organization of polymer brushes under poor solvent conditions. Polymer brushes in a poor solvent cannot undergo macroscopic phase separation due to the physical constraint imposed by the grafting points; therefore, they microphase separate forming aggregates. Under some conditions, the theory predicts that the homogeneous brush and the aggregates can exist as two different minima of the free energy. The theoretical methodology introduced in this work allows us to predict the minimum free-energy path connecting these two minima as well as the morphology of the system along the path. It is shown that the transition between the homogeneous brush and the aggregates may involve a free-energy barrier or be barrierless depending on the relative stability of the two morphologies and the chain length and grafting density of the polymer. In the case where a free-energy barrier exists, one of the morphologies is a metastable structure and, therefore, the properties of the brush as the quality of the solvent is cycled are expected to display hysteresis. The theory is also applied to study the adhesion/deadhesion transition between two opposing surfaces modified by identical polymer brushes and it is shown that this process may also require surpassing a free-energy barrier.

which is eq (13) in the main text.The right hand side of eq (S1) can be obtained from eq (9) in the main text: where   /, F P j   is given by eq (6) in the main text: ' The evaluation of   /, F P j   is the most expensive part of the calculation and therefore we optimized its numerical computation.We will explain next our approach for calculating eq.(S3).For simplicity, we will write the equation for a fixed j and, therefore, we will drop this index from our notation.In practice, the evaluation of   /, F P j   needs to be performed for all j, but this process can be trivially parallelized.In our code, the parallelization in j is distributed among processors using the MPI (Message Passing Interface) standard.
We first discretize our equations into a lattice.Since our problem is two-dimensional, we discretize our system in the x and y directions in square cells of side  (where  is a discretization parameter that we choose to be equal to 0.5 nm in the present work).The discretized system has size Dx in the x dimension and Dy in the y dimension.While in principle two indexes can be used to label each cell (ix and iy, representing the coordinates of the cell in x and y dimensions), the matrix formulation of the problem requires to use only one index to label each cell.Thus, we define this index as i = ix + (iy-1)*Dx (note that 1 ≤ i ≤ DxDy).We now discretize eq (S3) as: where np(,i) is the number of polymer segments that the chain in conformation  has in the lattice cell i and g(i,i') is a coefficient that the gives the strength of the VdW interactions between a segment in cell i and other segment in cell i' (see next section for its calculation).
It is also useful to define:

S3
In order to formulate the problem in matricial form, we define the following vectors: and the matrix: We can now rewrite eq (S4) as: At this point it is important to note that Nconf = 2.510 5 and DxDy  400.Therefore the matrix np has  10 8 elements and the matrix multiplication in eq (S10) becomes one of the bottlenecks of the MFEP calculation.The calculation of pot(i) from eq (S5) is much faster than this matrix multiplication, but it requires to know p(i) and s(i) (the density of the polymer and solvent at cell i, respectively).The evaluation of p(i) is a time consuming process by itself as it also requires a matrix multiplication (s(i) can be straightforwardly determined once p(i) is known using the packing constraint, eq (3) in the main text).
To write down the calculation of p(i) in matricial form, we start with eq (2) in the main text: We can rewrite this equation as The calculation of p(r,j) can be therefore trivially parallelized, so we will drop the index j in eq.(S13) for simplicity.We discretized the resulting equation, which results in: We now define the vectors: (1) ... ()   (1) ,..., ( ) so that eq (S14) can be rewritten in matricial form as: where the superindex T denotes the transverse matrix.
The most expensive computation steps in the calculation of the MFEP are the matrix multiplications in eqs (S10) and (S17), which involve the matrix np.This matrix contains the number of segments that each polymer conformation has in each lattice site.For each chain in the system, we require to know the position of NconfL segments (where L is the chain length of the polymer), while np has a size equal to DxDyNconf.The maximum fraction of non-zero elements in np is, therefore, L/(DxDy).Using typical values of L  50 and Dx  Dy  20, we get a maximum fraction of non-zero elements of  0.125 and, thus, np is a sparse matrix (the real fraction of non-zeros is actually smaller than this because some conformations have more than one segment in a given lattice cell).We store np as a compressed matrix in memory due to its large size.In previous works, 1, 2 we compressed np using a matrix of size NconfL where the element at the matrix position (,k) contained the position in the lattice of the segment k of conformation .In this previous works, we performed the matrix multiplications required in eqs (S10) and (S17) using loops over L and Nconf.Note that this procedure uses indirection to perform the matrix multiplication.In the present work, we decided to optimize the matrix multiplication step by using the subroutine dcsrmv from the Intel® MKL library, 3 which allows highly efficient multiplication of sparse matrices compressed in the CSR column format.The MKL implementation of our present program runs between 30% and 260% faster (depending on the parameters of the calculation and the computational architecture) than an implementation using the previous method for matrix multiplication.

Determination of the Poor Solvent Coefficients, g(i,i')
The function g(|r -r'|) determines the distance dependence of the VdW interactions in the theory.We model this dependence with the attractive branch of the Lennard-Jones potential with a cutoff of 2.5 , namely: where a is the segment length (we used a = 0.5 nm).The discretization of g(|r -r'|) into the lattice results in the coefficients g(i,i'), which are required to calculate eq (S5).This coefficients are obtained by integration of g(|r -r'|): if we denote the x and y coordinates of cells i and i' in the two dimensional lattice as (ix, iy) and (ix', iy'), respectively, then g(i,i') is given by: , ' , ,

Numerical Confirmation that the Barriers in the MFEP are Saddle Points of the Free-Energy Functional
In this section, we propose a method to confirm that the free-energy barriers of the MFEP found by the string method/molecular theory are actually saddle points of the free energy of the system.We have shown that taking This equation is valid not only for the local minima of the free energy, but also for all other stationary points (saddle points); although in practice numerical solvers fail to find the saddle points (which are unstable) and provide always the local minima of F (the stable solutions).We can now replace eq.(S21) into the free-energy functional of the system, eq (5) in the main text.Combining the resulting expression with the definition of the polymer number density (eq (2) in the main text), the packing constraint (eq (3) in the main text), the probability normalization constraint (eq (4) in the main text) and the definition of (r) (eq (7) in the main text) results in: 4 where the superindex "stationary" denotes that this expression for the free energy is only valid in the stationary points of F, in contrast with the general expression given by eqs (1) and (5) in the main text that is always valid.Therefore, if we now plot F and F stationary along the MFEP, we should find that F stationary = F only for the points that are S6 saddle points or extrema of F. Figure S1 shows that this is actually what happens for the initial and final states (that are local free-energy minima) and for the barrier state of the MFEP which is, therefore, a saddle point of F.
Figure S1: A. Free energy of the system, F, determined with eq (5) in the main text and free energy for the stationary points, F stationary , given by eq (S22) as a function of the normalized arc-length of the string (same data as in Figure 3 of the main text).B. Difference between the free energies in A as a function of the normalized arc-length of the string.