PSCF v1.1
Appendix: Self-Consistent Field Theory

Python Tools (Prev)         Appendix: Periodic Functions and Fourier Series (Next)

Below, we describe the self-consistent field theory (SCFT) equations that are solved by PSCF, giving a complete overview of the mathematical system of equations that are solved by the software and summarizing their connection to corresponding physical properties. This is not intended to be an introduction to the field theory itself, and we do not provide a derivation of the SCFT equations. We assume readers are familiar with the principles of polymer physics and statistical mechanics.

For readers looking to better understand the underlying theory and its applications to polymer science, we recommend the following references:

  • Fredrickson, Glenn H. The Equilibrium Theory of Inhomogeneous Polymers, Oxford University Press, 2006.
  • Arora, A., et al. "Broadly Accessible Self-Consistent Field Theory for Block Polymer Materials Discovery," Macromolecules, 2016, 49, 4675.
  • Matsen, M. W. "The standard Gaussian model for block copolymer melts," J. Phys.: Condens. Matter, 2001, 14, R21.
  • Shi, A.-C. "Self-Consistent Field Theory of Inhomogeneous Polymeric Systems: A Variational Derivation," Adv. Theory Simul., 2019, 2, 1800188.

Model and Notation

PSCF solves the SCFT equations for an incompressible mixture of any number of block polymer species and point-like solvent molecular species.

In what follows, we consider a system with \( C \) distinct monomer types, \( P \) block polymer species, \( S\) solvent species. We use a convention in which the symbols \( \alpha, \beta = 0, \ldots, C-1\) are used as indices for monomer types, symbols \( a, b \) are used to denote indices for blocks within a specific block polymer species, and symbols \( i, j\) are used to denote molecular species indices. In these notes, we use a convention in which species indices are ordered with all polymeric species listed first, so that species index values in the range \( i, j = 0, \ldots, P-1\) denote polymeric species, and values in the range \( P,\ldots, P+S-1 \) denote solvent species. A zero-based indexing scheme for indices is used throughout the notes to facilitate comparison to C/C++ source code.

To characterize inhomogeneous states, we consider the average number concentration \(\langle c_{\alpha}({\bf r}) \rangle\) for monomers of each type \( \alpha \) at each position \( {\bf r} \). In a model for an effectively incompressible liquid in which each monomer is assumed to occupy a volume \( v \), we may define a corresponding volume fraction

\[ \phi_{\alpha}(\textbf{r}) = v \langle c_{\alpha}({\bf r}) \rangle \quad, \quad\quad\quad\quad\quad (A.1) \]

for each monomer type, where \( v \) is a monomer reference volume that, by convention, has the same value for all monomer types. Incompressibility is imposed by requiring that

\[ 1 = \sum_{\alpha=0}^{C-1} \phi_{\alpha}({\bf r}) \quad\quad\quad\quad\quad (A.2) \]

at every position \( {\bf r} \) within the problem domain.

Self-Consistent Field Approximation

SCFT for a liquid of flexible polymers is based on a mean-field approximation that allows us to approximately predict properties of an interacting liquid by considering the behavior of a corresponding inhomogeneous gas of noninteracting molecules. We refer to this hypothetical non-interacting gas as the ideal-gas reference system. The ideal gas reference system must have the same chemical composition as the real system of interest and an inhomogeneous average concenration \( \langle c_{\alpha}({\bf r}) \rangle \) for each monomer type \( \alpha \) that is equal to that in the real system. Inhomogeneous monomer concentrations are induced in this reference system by a set of potential energy fields that couple to different types of monomer: Each monomer of type \( \alpha \) at position \( {\bf r} \) within this ideal gas is subjected to a potential \( k_{B}Tw_{\alpha}(\textbf{r}) \). In a self-consistent field approximation, this potential represents a free energy cost for placing a monomer of a specified type at a specified position, due to interactions among monomers within an inhomogeneous structure.

In SCFT, the monomer chemical potential fields (or "w fields") \( w_{0}, \ldots, w_{C-1} \) are assumed to depend upon volume fraction fields in a manner that reflects the effects of interaction among species in the mixture. PSCF currently uses a simple approximation for this dependence that is based on an inhomogeneous generalization of the parameterization of excess free energy used in the Flory-Huggins theory of homogeneous polymer mixtures. In this approximation, a solution to the SCFT problem is required to satisfy the self-consistent field (SCF) equations

\[ w_{\alpha}({\bf r}) = \sum_{\beta=0}^{C-1} \chi_{\alpha\beta} \phi_{\beta}({\bf r}) + \xi({\bf r}) \quad\quad\quad\quad\quad (A.3) \]

for all \( \alpha = 0,\ldots,C-1 \), in which \( \chi_{\alpha\beta}\) is a binary Flory-Huggins dimensionless interaction parameter, defined such that \( \chi_{\alpha\beta} = \chi_{\beta\alpha} \). Here, \( \xi({\bf r}) \) is a dimensionless Lagrange multiplier pressure field that must be chosen so that the monomer concentrations in the ideal gas reference system satisfy the incompressibility constraint given above.

In SCFT, average monomer concentration and volume fraction fields are computed by computing these properties for the ideal gas reference system. An iterative solution of an SCFT problem thus starts with an initial guess for the w fields for all monomer types, followed by iterative adjustment of these fields until the w fields and average monomer concentrations (or volume fractions) computed using these fields as inputs to the ideal gas model satisfy the above SCF equations and the incompressibility constraint.

Polymer Species

Each polymer species in a system modeled by PSCF must be a block polymer consisting of one or more blocks connected at junction points. Homopolymers are treated as a special case of a block polymer with only one block. Let \( B_{i} \) denote the number of blocks in block polymer species number \( i \), for \( i = 0, \ldots, P-1 \).

Branched Polymer Topology

Each block within a block polymer connects two vertices. Each such vertex is either a junction to which two or more blocks are connected or a chain end to which only one block is attached. In the language of graph theory, every allowed polymer structure corresponds to an acyclic connected graph (also known as a "tree") with some number of "edges" (or blocks) and "nodes" (or vertices). PSCF can only handle polymers with the topology of a connected acyclic graph, in which there is exactly one path (i.e., sequence of connected edges) connecting any two vertices. A standard theorem of graph theory aserts that any connected acyclic graph with \( B_{i} \) edges must have exactly \( B_{i} + 1 \) vertices or nodes. An acyclic block polymer with \( B_{i} \) blocks thus always has exactly \(B_{i} + 1\) vertices (i.e., junctions and free chain ends). For example, a linear triblock terpolymer with three blocks in a row has four vertices (two free ends and two junctions), and so does a star triblock copolymer with three blocks that meet at a central junction (three free ends and one junction).

PSCF assigns each block within a polymer species \(i\) with \( B_{i} \) blocks an integer index in the range \( [0, \ldots, B_{i}-1 ]\), and assigns each vertex (each junction or chain end) an index value in the range \( [ 0, \ldots, B_{i} ]\). The structure of a polymer can thus be characterized by assigning each block a monomer type, a contour length, and a pair of integer ids for the two vertices to which the block is attached. The length of each block is given by the ratio of the volume occupied by the block to the value of a monomer reference volume. Let \( L_{ia} \) denote the length of block \( a \) of polymer species \( i \).

Propagator Functions

For each block, we define two "propagator" functions, denoted generically by \( q( \textbf{r}, s) \), each of which is associated with one of the two directions in which the block may be traversed. Each propagator is said to originate at one of the two vertices to which the associated block is attached and terminate at the other. The contour variable parameter \( s \) is defined for each propagator \( q(\textbf{r},s) \) such that \( s = 0 \) at the vertex at which the propagator originates and \( s = L \) at the vertex at which the propagator terminates, where \( L \) is the contour length of the block. When we need to distinguish between the propagators associated with a single block, we will refer to them instead as \( q \) and \( q^{\dagger} \), though the choice of which is referred to as \( q \) and which \( q^{\dagger} \) is arbitrary.

For each vertex that is connected to \( M \) blocks, there are \( M \) propagators that originate at that vertex, which we refer to the "outgoing" propagators for that vertex, and \( M \) propagators that terminate at that vertex, which we refer to as "incoming" propagators. In what follows, we refer to a vertex with \( M \geq 2 \) attached blocks as an M-fold junction, and a vertex with only one attached block as a chain end.

It is useful for some purposes to imagine a process in which we divide a block of length \( L \) into two blocks of length \( s \) and \( L - s \) by adding a "virtual" two-fold junction at some coordinate \( s \) within the block. This yields a two-fold vertex that is connected to two blocks of the same monomer type. The insertion or removal of this type of virtual vertex modifies the graph that we use to describe a polymer but does change our physical description of the underlying underlying molecule.

Modified Diffusion Equation (MDE)

Each propagator \( q(\textbf{r}, s) \) associated with a block of monomer type \( \alpha \) can be shown to obey a partial differential equation (PDE)

\[ \frac{\partial q}{\partial s} = -H_{\alpha} q \quad\quad\quad\quad\quad (A.4) \]

in which \( H_{\alpha} \) is a linear differential operator whose action on any function \( f({\bf r}) \) of position is defined to be

\[ H_{\alpha} f({\bf r}) = - \frac{b_{\alpha}^{2}}{6} \nabla^{2} f({\bf r}) + w_{\alpha}(\textbf{r}) f({\bf r}) \quad\quad\quad\quad\quad (A.5) \]

in which \( b_{\alpha} \) is the statistical segment length for monomers of type \( \alpha \) , and \( w_{\alpha}\) is the dimensionless SCF chemical potential field for such monomers. We refer to this PDE as the modified diffusion equation (MDE). Numerical solution of the MDE is the most computationally expensive operation in any SCFT program.

Initial Conditions :

The initial condition for propagators that originate at a chain end requires that

\[ q(\textbf{r},s=0) = 1 \quad\quad\quad\quad\quad (A.6) \]

for all \( \textbf{r} \) in the relevant spatial domain.

The initial value \( q(\textbf{r}, s=0) \) of an outgoing propagator that originates at an \( M \)-fold junction is given by a pointwise product of values of the \( M - 1 \) incoming propagators for that junction, excluding the incoming propagator associated with the same block as that of the outgoing propagator of interest. Suppose we denote the outgoing propagators associated with the M blocks attached to an M-fold junction \( q_{a}({\bf r}, s) \), and denote the final or "tail" values of the M incoming propagators by \( q^{\dagger}_{a}({\bf r}) \), where \( a = 0, \ldots, M-1 \) is an integer identifer for the attached blocks. The initial condition for any propagator that originates at this vertex is then

\[ q_{a}({\bf r}, 0) = \prod_{b \neq a} q^{\dagger}_{b}({\bf r}) \quad\quad\quad\quad\quad (A.7) \]

where the product over \( b \) is taken over all \([0, \ldots, M-1] \) except the value \( a \) corresponding to the block associated with the outgoing propagator of interest.

In the simple case of a linear polymer, for which every junction connects two blocks, the initial condition for the each outgoing propagator that originates at a junction is simply equal to the final value of the incoming propagator associated with the other block that is attached to that junction. Consider a linear polymer with \( B \) blocks in which blocks and vertices are indexed sequentially by their position relative to one end of the chain. This is the indexing scheme for blocks and vertices that is used internally by PSCF for polymers that are assigned a linear type parameter value within the parameter file, which is the default parameter value. Let \( L_{a} \) denote the length of block \( a \), which connects vertices \( a \) and \( a + 1 \). Let \( q_{a} \) denotes the propagator in block \( a \) that originates at vertex \( a \), and let \( q_{a}^{\dagger} \) denotes the propagator in block \( a \) that originates at vertex \( a+1 \). With this notation, \( q_{a+1}(\textbf{r},0) = q_{a}(\textbf{r},L_{a}) \) and \( q_{a}^{\dagger}(\textbf{r},0)= q_{a+1}^{\dagger}(\textbf{r},L_{a+1})\) for \( a = 0, \ldots, B - 2 \).

As an example of a branched polymer, consider the case of a three armed star containing blocks with indices 0, 1 and 2 that are all connected at a central junction. Let \( L_{a} \) denote the length of block \( a \). Let \( q_{a}(\textbf{r},s)\) denote the propagator in block \( a \) that originates at the central junction, and let \( q_{a}^{\dagger}(\textbf{r},s) \) denote the propagator in block \( a \) that originates at a free end and terminates at the central junction. The initial condition for propagators that originate at chain ends requires that \( q_{a}^{\dagger}(\textbf{r},s=0) = 1 \) for all \( \textbf{r} \) and all \( a \). The initial condition for the outgoing propagator in block 0 that originates at the central vertex is then given by, for example, by

\[ q_{0}(\textbf{r}, s=0) = q_{1}^{\dagger}(\textbf{r}, L_1) q_{2}^{\dagger}(\textbf{r}, L_2) \quad. \quad\quad\quad\quad\quad (A.8) \]

Initial conditions for \( q_{1} \) and \( q_{2} \) are analogous.

Boundary Conditions :

In the pscf_pc and pscf_pg programs, the MDE for propagators in a periodic structure is solved within a single unit cell, subject to periodic boundary conditions.

In the pscf_fd program, the MDE for each propagator \( q(x,s) \) is solved within a range of values of the single spatial coordinate \( x \), where \( x \) may represent either a Cartesian normal coordinate for a system with the geometry of a flat film or a radial coordinate for a system with cylindrical or spherical symmetry. In the case of a problem that is defined in cylindrical or spherical coordinates with a domain that includes \( x = 0 \), we must require that \( \partial q/\partial x = 0 \) at \( x = 0 \) in order to avoid introducing gradient discontinuities at the origin. In most foreseeable applications of this program, one or both of the limits of the computational domain is intended to represent a point that lies deep within homogeneous region, such as bulk region surrounding a micelle or a bulk phase far from an interface. In these cases, it is also appropriate to set \( \partial q/\partial x = 0 \) at these boundaries. In light of these considerations, the pscf_fd program is designed to apply von Neumann boundary conditions requiring that \( \partial q/\partial x = 0 \) at both upper and lower limiting values of \( x \).

Monomer Probability Densities (Unnormalized)

Values of propagator functions that are obtained by solving the modified diffusion equation can be used to compute the probability of finding any point on a polymer chain within an ideal gas reference system at any at specified position within the relevant spatial domain: The probability per unit volume of finding a specified vertex of a polymer at a specific local \( {\bf r} \) can be shown to be proportional to the product of the incoming propagators for that vertex, at the the position \( {\bf r} \).

This relationship applies to an optional or "virtual" 2-fold vertex inserted at any point within any block of a an acyclic branched polymer, as well as to vertices that are necessary to describe the polymer topology. Let \( \rho({\bf r}, s) \) denote the probability density per unit volume that a monomer or point with a contour coordinate \( s \) within a specified block of length \( L \) will be found within an infinitesimal volume centered at a position \( {\bf r} \). This quantity is given to within a constant of proportionality by

\[ \rho({\bf r}, s) \propto q({\bf r}, s) q^{\dagger}({\bf r}, L-s) \quad, \quad\quad\quad\quad\quad (A.9) \]

Here, \( q({\bf r}, s) \) and \( q^{\dagger}({\bf r}, L-s) \) are the two propagators associated with that block, which are also the two incoming propagators associated with a virtual vertex introduced at that point within the block. The two propagators in this product are evaluated at contour coordinates that correspond to the same point within the block, but are each evaluated at a value \( s \) or \( L - s \) that corresponds to the contour distance measured from the vertex at which that propagator originates. The more general expression for the probability density for the position of an M-fold vertex involves the product of the M incoming propagators, all evaluated at the same position \( r \), and at contour coordinate values that correspond to the lengths of the associated blocks.

We show below that the constant of proportionality that relates the probability density for a vertex to the product of incoming propagators is the same for every vertex in a polymer, as well as for any possible virtual vertex (i.e., any monomer within any block of that polymer).

Theorems about Propagator Products

We may use the MDE and its associated initial and boundary conditions to prove the following series of theorems about products of incoming propagators associated with any vertex in a polymer.

Theorem A.1 : Suppose \( q({\bf r}, s) q^{\dagger}({\bf r},L-s)\) is the product of two propagators defined within a single block of length \( L \). The value of this product at \( s = 0 \), corresponding to its value at the vertex at which \( q \) originates, is equal to the product of all the incoming propagators for this vertex. Similarly, the value of this product at \( s = L \), corresponding to its value at the vertex at which \( q^{\dagger} \) originates, is equal to the product of all incoming propagators at this vertex.

Proof : The value \( q^{\dagger}(\textbf{r},L-s) \) at \( s = 0 \) yields the final value for the incoming propagator for the vertex at \( s = 0 \) that is associated with the block within which \( q \) and \( q^{\dagger} \) are both defined. The initial condition for \( q \) at \( s = 0 \) for an M-fold junction with \( M > 1 \) requires that \( q({\bf r}, s=0) \) is equal to the product of the final values of incoming propagators associated with all other blocks that are attached to that vertex , excluding the incoming propagator associated with the block along which \( q \) and \( q^{\dagger} \) are defined. For any such junction, the value of the product \( q({\bf r},s) q^{\dagger}({\bf r},L-s) \) at \( s = 0 \) is thus equal to the product of all incoming propagators at that vertex. If \( q \) originates at a chain end, then \( q(\textbf{r}, s) = 1 \), and so this product simply yields the value of the single incoming propagagator \( q^{\dagger}(\textbf{r}, L) \) at that chain end, which in this case is equal to the product of incoming propagators. The corresponding statement about the value of \( q q^{\dagger} \) at \( s = L \) may be obtained by completely analogous reasoning, with the roles of \( q \) and \( q^{\dagger} \) reversed.

Theorem A.2 : Let \( I(s) \) represent the spatial integral

\[ I(s) \equiv \int d\textbf{r} \; q (\textbf{r},s) q^{\dagger}(\textbf{r},L-s) \quad\quad\quad\quad\quad (A.10) \]

of a product of two propagators that are both defined within a block of length \( L \), where \( s \) is the contour variable for a point within that block. If the propagators \( q \) and \( q^{\dagger} \) obey a modified diffusion equation within some domain, subject to periodic, Dirichlet, von Neumann or Robin boundary conditions, then the value of \(I(s)\) is independent of contour coordinate \( s \) for all \( s \in [0,L] \).

Proof : Note that

\begin{eqnarray*} \frac{dI(s)}{ds} & = & \frac{d}{ds} \int d\textbf{r} \; q (\textbf{r},s) q^{\dagger}(\textbf{r},L-s) \\ & = & \int d\textbf{r} \; \left \{ \frac{\partial q (\textbf{r},s)}{\partial s} q^{\dagger}(\textbf{r},L-s) - q (\textbf{r},s) \frac{\partial q^{\dagger}(\textbf{r},L-s)}{\partial s} \right \} \\ & = & -\int d\textbf{r} \; \left \{ [ \hat{H} q (\textbf{r},s) ] q^{\dagger}(\textbf{r},L-s) - q (\textbf{r},s) [\hat{H} q^{\dagger}(\textbf{r},L-s)] \right \} \quad. \quad\quad\quad\quad\quad (A.11) \end{eqnarray*}

We may then use the fact that the operator \( \hat{H} \) (like the analogous quantum mechanical Hamiltonian operator for a particle in a potential) can be shown to be a self-adjoint differential operator for any of the boundary conditions given in the theorem statement to show that

\[ \int d\textbf{r} \; [\hat{H} q (\textbf{r},s) ] q^{\dagger}(\textbf{r},L-s) = \int d\textbf{r} \; q (\textbf{r},s) [\hat{H} q^{\dagger}(\textbf{r},L-s) \quad, \quad\quad\quad\quad\quad (A.12) \]

implying a cancellation of the two terms in the integrand of the last expression for \( dI(s)/ds \). This thus gives

\[ dI(s)/ds = 0 \quad\quad\quad\quad\quad (A.13) \]

for all \( s \in [0,L] \), implying that \( I(s) \) is independent of \( s \) within each block.

Theorem A.3 : The value of the spatial integral of the product of incoming propagators associated with a vertex of an acyclic branched polymer has the same value for every vertex in the polymer, including any virtual 2-fold vertex inserted into any block of the polymer.

Proof : Theorem A.2 shows that the value of the spatial integral \( I(s) \) of the product of incoming propagators associated with a virtual 2-fold vertex inserted into any block of a polymer is independent of the coordinate \( s \) within the block. Combining this with theorem 1, we may show that the integral of the product of incoming propagators must have the same value for any vertices that are directly connected by a block. Because all vertices in a polymer are connected by some sequence of blocks, this implies that the value of this spatial integral must be the same for every vertex in any acyclic branched polymer, including virtual vertices that we insert at any point within any block.

Molecular Partition Function

Let \( Q_{i} \) denote a normalized single-chain partition function for a polymer molecule of species \( i \) within an ideal gas reference system that is subjected to a specified set of monomer potential fields. We define the value of \( Q_{i} \) to be equal to the spatial average (i.e., the integral over the entire spatial domain of the problem divided by the volume of this domain) of the product of the values of the incoming propagators evaluated at any vertex of the polymer. The statement of theorem A.3 implies that we should get the same value for \( Q_{i} \) by evaluating it using any vertex in a polymer of species \( i \). Different values for \( Q_{i} \) are, however, generally obtained for different species within a mixture.

A particularly compact expression for this quantity may be obtained by evaluating the incoming propagator at any chain end, for which the vertex is attached to a single block of length \( L \). In this case, the above definition yields

\[ Q_{i} = \frac{1}{V} \int \! d\textbf{r} \; q^{\dagger}(\textbf{r},L) \quad\quad\quad\quad\quad (A.14) \]

where \( q^{\dagger}(\textbf{r}, L) \) denotes the final value of the single incoming propagator at that chain end. Here, the integral is taken over the entire spatial domain for the problem of interest, while \( V \) is the generalized volume of that domain. In the case of a periodic structure, the required spatial average may be evaluated by evaluating the integral over one crystallographic unit cell and then dividing by the generalized volume of that unit cell (i.e., the volume per unit cell for a 3D crystal, area per 2D unit cell for a 2D structure such as a hexagonal cylinder phase, or the length per unit cell in a 1D lamellar phase).

Another potentially useful expression for \( Q_{i} \) may be obtained by considering a virtual 2-fold vertex inserted at contour coordinate \( s \) within any block \( a \) of length \( L_{ia} \) of a polymer species of \( i \) and using the result of theorem A.3 to show that

\[ Q_{i} = \frac{1}{V} \int d\textbf{r} \; q (\textbf{r},s) q^{\dagger}(\textbf{r},L-s) \quad\quad\quad\quad\quad (A.15) \]

where \( L \) is the length of the chosen block, and where theorem A.3 guarantees that the value of the RHS of this equation is actually independent of the choice of values for \( s \) and block index \( a \).

Monomer Probability Densities (Normalized)

We can use the above definition of \( Q_{i} \) to write an explicit expression for the absolute probability density for the location of any specified monomer within an ideal gas reference system. Let \( \rho_{i a}({\bf r}, s) \) denote the probability density for the position \( {\bf r} \) of a monomer with contour position \( s \) with block \( a \) of a single polymer of species \( i \) within a system that contains only a single polymer of this species, with no non-bonded pair interactions. By combining our previous proportionality that gives this quantity to within a constant with the requirement that integral of this probability density over the problem domain must equal unity, we find that

\[ \rho_{ia}({\bf r}, s) = \frac{1}{I_{ia}(s)} q({\bf r}, s) q^{\dagger}({\bf r}, L_{ia} - s) \quad, \quad\quad\quad\quad\quad (A.16) \]

in which \( q \) and \( q^{\dagger} \) are the two propagators associated with the block of interest, \( L_{ia} \) is the length of this block, and \( I_{ia}(s) \) is the spatial integral defined in the statement of theorem A.2, defined here for monomer with contour coordinate \( s \) within block \( a \) of polymer species \( i \). Using the expression given above to relate \(Q_{i}\) and \( I_{ia}(s) \), we find that

\[ \rho_{ia} ({\bf r}, s) = \frac{1}{V Q_{i}} q({\bf r}, s) q^{\dagger}({\bf r}, L-s) \quad\quad\quad\quad\quad (A.17) \]

for any choice of block index \( a \) and contour coordinate \( s \).

Block Volume Fractions

Consider a system that contains \( M_{i} \) chains of polymer species \( i \) within a total volume \( V \). The volume occupied by each such polymer is \( v L_{i} \), where \( L_{i} \) is the total length of polymer species \( i \), given by the sum of the length of all of the blocks of that polymer. Let \( \overline{\phi}_{i} \) denote the average overall volume fraction occupied by molecule species \( i \) within the mixture. This quantity is given by a ratio

\[ \overline{\phi}_{i} = \frac{M_{i} L_{i} v}{V} \quad\quad\quad\quad\quad (A.18) \]

where the numerator of the RHS is the total volume occupied by \( M_{i} \) such chains.

Let \( \phi_{ia}(\textbf{r}) \) denote the local volume fraction of monomers belonging to block \( a \) of polymer species \( i \). This quantity is related to the probability density \( \rho_{ia}({\bf r},s) \) defined above by

\[ \phi_{ia}({\bf r}) = v M_{i} \int\limits_{0}^{L_{ia}} \! ds \; \rho_{ia}({\bf r}, s) \quad, \quad\quad\quad\quad\quad (A.19) \]

where \( M_{i} \) is the total number of chains of species \( i \) and \( L_{ia} \) is the length of the relevant block. Rewriting \( v M_{i} \) in terms of \( \overline{\phi}_{i} \) then yields

\[ \phi_{ia}(\textbf{r}) = \frac{\overline{\phi}_{i}}{L_{i}Q_{i}} \int\limits_{0}^{L_{ia}} \! ds \; q(\textbf{r},s) q^{\dagger}(\textbf{r},s) \quad\quad\quad\quad\quad (A.20) \]

where \( L_{i} \) is the total contour length of all of the blocks of polymer species \( i \).

The total volume fraction \( \phi_{\alpha}(\textbf{r}) \) for each monomer type \( \alpha\) is simply the sum of volume fractions for all polymer blocks and all solvent species that have been assigned a monomer type index \( \alpha\) .

Solvent Species

Each solvent species \( i\) is associated with a specific monomer type \( \alpha\) and a dimensionless "size" parameter. The size parameter for a solvent is a measure of volume closely analogous to the length parameter for a block in a block polymer, and is equal to the ratio of the volume per molecule occupied by the solvent to the monomer reference volume. In light of this analogy, we denote the "size" of solvent \( i\) here by \( L_{i} \).

In the single molecule problem for solvent species, the free energy penalty for a solvent molecule of monomer type \( \alpha\) to be located at position \( \textbf{r} \) is given by \( k_{B}T L_{i} w_{\alpha}(\textbf{r}) \) . This implies the existence of a solvent volume fraction field \( \phi_{i}(\textbf{r}) \propto \exp[-L_{i}w_{\alpha}(\textbf{r})] \) .

We define a normalized molecular partition function \( Q_{i} \) for solvent species \( i \) that is given by an integral

\[ Q_{i} = \frac{1}{V} \int \! d\textbf{r} \; \exp (-L_{i}w_{\alpha}(\textbf{r}) ) \quad\quad\quad\quad\quad (A.21) \]

for all \( i = P, \ldots, P + S - 1\) . The local volume fraction \( \phi_{i}(\textbf{r}) \) for solvent species \( i \) is then given by a ratio

\[ \phi_{i}(\textbf{r}) = \frac{\overline{\phi}_{i}}{Q_{i}} \exp(-L_{i}w_{\alpha}(\textbf{r})) \quad\quad\quad\quad\quad (A.22) \]

in which \( \overline{\phi}_{i} \) is the overall volume fraction of solvent species \( i \) within the mixture. Note that the prefactor in this expression is chosen to guarantee that the spatial average of \( \phi_{i}(\textbf{r})\) (the integral divided by total volume) is give by the parameter \( \overline{\phi}_{i} \).

Open and Closed Ensembles

PSCF allows the user to choose either an open or closed statistical ensemble for each molecular species in the mixture, by allowing a user to enter a value for either the overall volume fraction for that species, denoted here by \( \overline{\phi}_{i} \), or a value for the chemical potential, denoted here by \( \mu_{i} \).

To use an open ensemble for a particular species in SCFT, one uses an open (i.e., grand-canonical) ensemble to treat the behavior of that species when considering the hypothetical ideal gas of non-interacting molecules subjected to a potential. Use of a grand-canonical ensemble to describe to a gas of molecules of species \( i \) yields an overall concentration for that species and a corresponding overall volume fraction \( \overline{\phi}_{i} \) that are both proportional to the product \( e^{\mu_{i}/k_{B} T} Q_{i} \). The value of the constants of proportionality in the expressions depends upon arbitrary choices of convention for normalizing volume used to relate the configurational integral for a single chain to a dimensionless partition function and a choice of reference state used to define the chemical potential. In the equations used in PSCF, we set

\[ \overline{\phi}_{i} = e^{\mu_{i}/k_{B} T} Q_{i} \quad, \quad\quad\quad\quad\quad (A.23) \]

as a matter of convention.

Calculations of the local volume fractions for monomers within a specific block of a block copolymer and for the local volume fraction of a solvent species are thus performed for molecular species that are treated in an open ensemble, by simply replacing the parameter \( \overline{\phi}_{i} \) by the above expression in all expressions that involve \( \overline{\phi}_{i} \).

The value of the chemical potential is computed for any species that is treated in a closed ensemble (i.e., any species for which a volume fraction is specified in the parameter file) by computing

\[ \frac{\mu_{i}}{k_{B}T } = \ln \left ( \overline{\phi}_{i} / Q_{i} \right ) \quad, \quad\quad\quad\quad\quad (A.24) \]

after computing \( Q_{i} \).

Values of species chemical potentials in both input and output files used by PSCF, which are denoted indicated the label "mu", are given in thermal energy units, and thus are correspond to values of the ratio \( \mu_{i}/k_{B}T \).

Thermodynamic Properties

PSCF reports values of the Helmholtz free energy and the pressure in non-dimensionalized form.

The SCFT prediction for the total Helmholtz free energy \( F \) can be expressed as a sum

\[ F = F_{\rm id} + F_{\rm int} \quad\quad\quad\quad\quad (A.25) \]

where \( F_{\rm id} \) is a "ideal" contribution that arises from polymer configurational entropy, and \( F_{\rm int} \) is an "interaction" contribution that arises directly from excess free energy of interaction that is characterized here in terms of binary interaction parameters. The non-dimensionalized free energy that is reported as "fHelmholtz" in the PSCF thermo file block format is the free energy per monomer reference volume in thermal energy units, also denoted by the symbol

\[ \tilde{f} \equiv \frac{F v}{ V k_{B} T } \quad. \quad\quad\quad\quad\quad (A.26) \]

This quantity can be expressed as a sum

\[ \tilde{f} = \tilde{f}_{\rm id} + \tilde{f}_{\rm int} \quad\quad\quad\quad\quad (A.27) \]

in which the ideal and interaction contributions are given by

\begin{eqnarray*} \tilde{f}_{\rm id} & = & \sum_{i=0}^{P+S-1} \frac{\overline{\phi}_{i}}{L_{i}} \left [ \ln ( \overline{\phi}_{i} / Q_{i} ) - 1 \right ] - \frac{1}{V} \sum_{\alpha=0}^{C-1} \int \! d\textbf{r} \; w_{\alpha}(\textbf{r}) \phi_{\alpha}(\textbf{r}) \\ \tilde{f}_{\rm int} & = & \frac{1}{2V} \sum_{\alpha, \beta =0}^{C-1} \chi_{\alpha\beta} \int \! d\textbf{r} \; \phi_{\alpha}(\textbf{r}) \phi_{\beta}(\textbf{r}) \quad\quad\quad\quad\quad (A.28) \end{eqnarray*}

Note that the sum over species in the first term in the definition of \( F_{\rm id} \) is a sum over all species, including polymeric and solvent species. Unified notation for both types of molecule has been made possible by using the symbol \( L_{i} \) to denote a ratio of total molecular volume to monomer reference volume for species \( i \). This quantity is given by the overall length per molecule for block polymer species (the sum of the lengths of all blocks), or by the "size" parameter for solvent species.

The macroscopic physical pressure \( P \) is computed from the identity

\[ P = - \frac{F}{V} + \sum_{i=0}^{P+S-1} \frac{\mu_{i}\overline{\phi}_{i}}{L_{i}v} \quad\quad\quad\quad\quad (A.29) \]

in which \( F / V \) is the Helmholtz free energy density and \( \overline{\phi}_{i}/(L_{i}v) = M_{}/V\) is the average number of molecules of species \( i \) per unit volume. The value given in the thermo block output file format, with the label "pressure", is the dimensionless quantity

\[ \tilde{p} \equiv \frac{Pv}{k_{B}T} = - \tilde{f} + \frac{\mu_{i}\overline{\phi}_{i}}{L_{i}k_{B}T} \quad. \quad\quad\quad\quad\quad (A.30) \]

In Eqs. (A.29) and (A.30), the sum over species represents a sum over both polymeric and solvent species, and \( L_{i} \) denotes a ratio of molecular volume to monomer reference volume, exactly as in Eq. (A.28).

Indeterminancy of Pressure in Canonical Ensemble

In grand-canonical and mixed ensembles, values for the Lagrange multiplier field \( \xi(\textbf{r}) \) and the macroscopic pressure \( P\) are uniquely determined by the values input for the chemical potentials.

In the special case of canonical ensemble, however, the value of the Lagrange multplier field \( \xi(\textbf{r}) \) is not uniquely determined by the SCF equations and incompressibility constraint. In PSCF, a calculation is performed in canonical ensemble if the parameter file provides a volume fraction, rather than a chemical potential, as an input parameter for every polymer and solvent species. In this case, the relevant working equations can be shown to define \( \xi(\textbf{r}) \) only to within shifts by an arbitrary spatially homogeneous constant. Change in the value of \( \xi(\textbf{r}) \) by a homogeneous constant cause corresponding changes in values computed for \( w \) fields, chemical potentials and the non-dimensional pressure \( \tilde{p} \). As a result, values for these quantities are also not uniquely defined in this ensemble, unless an additional mathematical constraint is imposed by convention to obtain a unique solution for \( \xi({\bf r}) \).

Different conventions for this additional constraint on \( \xi({\bf r}) \) are used in different PSCF programs. The implementations of the pscf_pc and pscf_pg programs for periodic microstructures obtain equations with a unique solution by requiring, as a matter of convention, that the spatial average of \( \xi(\textbf{r}) \) must vanish in canonical ensemble. The pscf_fd code for one-dimensional problems instead uses a convention that requires that, in canonical ensemble, the value of \( w_{\alpha}({\bf r}) \) must equal zero for the last monomer type ( \( \alpha = C-1 \)) at the last grid point.

When an SCFT calculation is performed in canonical ensemble, each PSCF program outputs converted values for the pressure, chemical potentials, and \( w \) fields that are all consistent with the arbitrary convention used in that program to define a unique solution for \( \xi \). Numerical results for these quantities are thus sensitive to this choice of convention. Results for the Helmholtz free energy density and volume fraction fields can, however, be shown to be independent of changes in the value of \( \xi \) by a homogeneous constant, and are thus independent of this choice of convention.

Results of a calculation that is initially carried out in canonical ensemble can be exactly reproduced in a subsequent calculation that is performed in a mixed or grand-canonical ensemble by using the values of some or all of the chemical potentials that were output by the canonical calculation as inputs to a calculation performed in mixed or grand-canonical ensemble. Attempting to reverse procedure, however, leads to more complicated relationship. Using the converged species volume fractions output by a calculation that is performed in a mixed or grand-canonical as inputs to a subsequent calculation that is performed in canonical ensemble will generally lead to results for values for the \( w \) fields, species chemical potentials and pressure that are different from those reported in the original calculation, because the value of \( \xi(\textbf{r}) \) will be shifted by a homogeneous constant in order to satisfy the convention that is imposed only in the special case of a canonical ensemble. Equal values for monomer volume fraction fields and Helmholtz free energies should, however, be obtained in any ensemble.

Overview: Working Equations, Inputs and Outputs

The above discussion provide a complete set of working equations for SCFT of an incompressible liquid in the form implemented in PSCF, for calculations in a fixed spatial domain. We simplify the discussion by restricting ourselves here to calculation that are carried out in a fixed one-dimensional or periodic spatial domain. Extension to SCFT calculations in a periodic unit cell with flexible parameters, as performed by pscf_pc or pscf_pg programs, would require an additional set of working equations to impose a requirement that the free energy density be minimized with respect to variation of the unit cell parameters. The working equations for a PSCF SCFT calculation in a fixed spatial domain are:

  • The modified diffusion equation [Eqs. (A.4) and (A.5)]
  • A definition of \( Q_{i} \) for each species [Eq. (A.14) for polymers and Eq. (A.21) for solvents]
  • The relationship between \( \overline{\phi}_{i} \) and \( \mu_{i} \) [Eq. (A.23) or (A.24)]
  • Expressions for block and solvent monomer volume fractions [Eq. (A.20) for a polymer block and (A.22) for a solvent]
  • Self-consistent field condition [Eq. (A.3)]
  • Incompressibility constraint [Eq. (A.2)]
  • Expression for Helmholtz free energy per monomer \( \tilde{f} \) [Eqs. (A.26 - A.28)]
  • Expression for non-dimensionalized pressure \( \tilde{p} \) [Eq. (A.30)]

The key outputs from solution of this set of equations are:

  • A converged monomer chemical potential fields \( w_{\alpha}({\bf r}) \) for each monomer type \( \alpha = 0, \ldots, C-1 \)
  • A converged monomer volume fraction field \( \phi_{\alpha}({\bf r}) \) for each monomer type \( \alpha \)
  • A computed chemical potential \( \mu_{i} \) or overall volume fraction \( \overline{\phi}_{i} \) for each polymer and solvent species \( i \) (depending on whether \( \overline{\phi}_{i} \) or \( \mu_{i} \) was provided as an input parameter for that species)
  • Values for non-dimensionalized Helmholtz free energy \( \tilde{f} \) and pressure \( \tilde{p} \).

The required inputs to such an SCFT calculation are:

  • A statistical segment length for each monomer type
  • A description of each block polymer type as a list of block lengths and vertex indices for attached vertices
  • A size (volume relative to monomer reference volume) for each monomer type
  • Values for all elements of the symmetric matrix of Flory-Huggins interaction parameters
  • A description of the spatial domain and boundary conditions

It is straightforward to confirm, however, that the monomer reference volume \( v \) does not appear as a parameter in any of the final working equations enumerated above, and is thus not a required input parameter for an SCFT calculation.

PSCF allows the user to optionally input a value for \( v \) in the parameter file (labelled as "vMonomer") because a value for \( v \) can be used for some calculations other than the standard SCFT calculation described above. For example, a value of \( v \) is required as an input to compute the total excess free energy of a spherical micelle, which can be output by the pscf_fd program when used to simulate this geometry. A value of \( v \) will also be required as an input for stochastic field theoretic simulations that are currently under development as a planned addition to PSCF.


Python Tools (Prev)         User Guide (Up)         Appendix: Periodic Functions and Fourier Series (Next)