On the Ternary Ohta–Kawasaki Free Energy and Its One-dimensional Global Minimizers

We study the ternary Ohta–Kawasaki free energy that has been used to model triblock copolymer systems. Its one-dimensional global minimizers are conjectured to have cyclic patterns. However, some physical experiments and computer simulations found triblock copolymers forming noncyclic lamellar patterns. In this work, by comparing the free energies of the cyclic pattern and some noncyclic candidates, we show that the conjecture does not hold for some choices of parameters. Our results suggest that even in one dimension, the global minimizers may take on very different patterns in different parameter regimes. To unify the existing choices of the long-range coefficient matrix, we present a reformulation of the long-range term using a generalized charge interpretation and thereby propose conditions on the matrix in order for the global minimizers to reproduce physically relevant nanostructures of block copolymers.

All data generated or analyzed during this study are included in the supplementary material.

The authors would like to thank Professors Chong Wang, Xiaofeng Ren, An-Chang Shi, Juncheng Wei and Yanxiang Zhao for helpful discussions.

Authors and Affiliations

1. Department of Applied Physics and Applied Mathematics, Columbia University, New York, NY, 10027, USA

   Zirui Xu

2. Department of Applied Physics and Applied Mathematics, and Data Science Institute, Columbia University, New York, NY, 10027, USA

   Qiang Du

Corresponding author

Correspondence to Zirui Xu.

Communicated by Rustum Choksi.

Publisher's Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

This work is supported in part by the National Science Foundation DMS-2012562 and DMS-1937254.

A Underlying Mechanism of Interactions Between Generalized Charges

The decomposition (5.10) of \([f_{ij}]\) suggests a possible way to interpret the interactions between the generalized charges:

• Each charge consists of two sub-charges, as shown in Fig. 9. For example, a charge of type 1 consists of a sub-charge of type and a sub-charge of type .

• For \(i\ne j\), the potential energy due to the interaction between the pair of sub-charges at \(\vec x\) and at \(\vec y\) is \(f_{ij}\,G(\vec x,\vec y)\), and that between at \(\vec x\) and at \(\vec y\) is \(-f_{ij}\,G(\vec x,\vec y)\). There is no interaction within other pairs.

Decomposition of charges into sub-charges. , and balls represent charges of types 1, 2 and 3, respectively. Each charge consists of two (out of six) types of sub-charges (Color figure online)

Full size image

We further assume that there are van der Waals forces between charges, that the cohesive forces between charges of the same type are stronger than the adhesive forces between charges of different types, so that the charges behave like immiscible fluids in the thermodynamic limit at certain temperature, with the interfacial tensions being \(\{c_{ij}\}\). (For a general account of the interfacial tension, see, e.g., Birdi (2015).) To ensure the overall charge neutrality, there must be the same number of charges of each type (Fig. 9 shows four charges of each type). In accordance with the volume constraints, the volume ratio of charges of types 1, 2 and 3 is \(\omega _1:\omega _2:\omega _3\). On the continuum level, such a discrete particle system can be described by the ternary O–K free energy and thus serves as an intuitive analogy. This analogy naturally raises a question on how to arrange balls in 1-D in order to minimize the potential energy between charges, a question further discussed in Sect. 6.2.1 and Appendix B. The decomposition into simple interactions between sub-charges helps us answer this question for \(f_{12},f_{13},f_{23}\leqslant 0\) (see Proposition 7.6).

B Optimal Arrangement of Charged Balls in 1-D

1.1 B.1 Binary case

Given a positive integer n, consider a 1-D periodic cell [0, 1] packed with n balls of type \(A\) and n balls of type \(B\) , with unit amounts of positive and negative point charges at their centers, respectively. We assume that all the balls have the same radius \(\frac{1}{4n}\), and that their centers are located at \(\frac{k}{2n}\) for \(k = 1,2,\cdots ,2n\). Let \(u:\big \{\frac{k}{2n}\big \}_{k=1}^{2n}\rightarrow \pm 1\) represent the arrangement of the balls, with 1 and \(-1\) denoting \(A\) and \(B\) , respectively, then the total potential energy between charges can be written as

where G is given by (7.4). By Proposition 7.1, the alternating arrangement \(A\)\(B\)\(\cdots \;\) \(A\)\(B\) minimizes U, as shown in Fig. 10 (Throughout Appendix B,“minimize” refers to “globally minimize”.). This is a long-range variant of the 1-D antiferromagnetic Ising model without external fields, subject to the zero overall spin constraint.

An optimal arrangement for \(n=3\). and represent \(A\) and \(B\) balls with positive and negative point charges at their centers, respectively (Color figure online)

Full size image

Proposition 7.1

The minimizer of (7.1) is (up to translation) \(u\big (\frac{k}{2n}\big )=(-1)^{k-1}\) for \(k = 1,2,\cdots ,2n\).

Proof

(Inspired by Ren and Wei (2000, Proof of Proposition 3.1) Within an optimal arrangement, let us prove that l must be 1 for any segment like the following

By translational invariance, we can assume that the above segment occupies the first \(l\!+\!2\) sites \(\big \{\frac{k}{2n}\big \}_{k=1}^{l+2}\). By assumption, U should not decrease if we swap the \(A\) and \(B\) at the first two sites. Let \(u^*\) represent this optimal arrangement, then we have

Define the electrostatic potential

Then from (7.2), we know

From (7.4) and \(\sum _{m=1}^{2n}u\big (\frac{m}{2n}\big )=0\), we can see that the potential V is piecewise quadratic in x with the coefficient of \(x^2\) being 0 and thus is linear on every subinterval \([\frac{k-1}{2n},\frac{k}{2n}]\) for \(k=1,\cdots ,2n\). Define the electrostatic field \(E(x;u)=\mathrm {d}V(x;u)/\mathrm {d}x\), then E is piecewise constant in x and

Analogously we have

We also know for \(k=1,\cdots ,2n\!-\!1\),

Consequently, we have \(-\,l>(-1)-1\), that is, \(l=1\). \(\square \)

Remark 7.2

The above proof can be generalized to the case where different types of balls have different sizes, that is, the positions of the point charges are no longer uniform (i.e., \(x_k=\frac{k}{2n}\) for \(k=1,2,\cdots ,2n\)). Instead, we have

where \(x_0=0\) is identified with \(x_{2n}=1\) so that \(u(x_0)=u(x_{2n})\), and the radii of \(A\) and \(B\) balls are \(\frac{1+\omega }{4n}\) and \(\frac{1-\omega }{4n}\), respectively, for some \(\omega \in (-1,1)\backslash \{0\}\).

1.2 B.2 Ternary Case

Given a positive integer n, consider a 1-D periodic cell [0, 1] packed with n balls of type \(A\) , n balls of type \(B\) , and n balls of type \(C\) , labeled by 1, 2 and 3, respectively. Balls of type i are assumed to have the radius \(\omega _i/(2n)\). We also assume \(f_{ij}\,G(x,y)\) to be the potential energy between a ball of type i centered at x and a ball of type j centered at y. The total potential energy U is the sum of all the pairwise interactions:

where G is given by (7.4). The k-th ball, which is of type \(i_k\), is centered at \(x_k\), and the m-th ball, which is of type \(j_m\), is centered at \(y_m\). Therefore, we have the following relation

for \(k=1,2,\cdots ,3n\), with \(x_0=y_0=0\) and \(i_0=j_0=i_{3n}=j_{3n}\).

Conjecture 7.3

For any admissible interaction strength matrix \([f_{ij}]\) and any positive \(\{\omega _i\}\) satisfying \(\sum _i\omega _i=1\), the arrangement \(A\)\(B\)\(C\)\(\cdots \)\(A\)\(B\)\(C\) (i.e., \(i_k\equiv k\mod 3\)) minimizes (7.3), as illustrated on the right side of Fig. 5.

We have numerically verified Conjecture 7.3 from \(n=2\) to 8 for the Cartesian product of 100 choices of \([f_{ij}]\) and 72 choices of \(\{\omega _i\}\) (up to permutations we can assume \(\omega _1\leqslant \omega _2\leqslant \omega _3\)), as shown in Figs.11 and 12. We also prove some special cases of Conjecture 7.3 in Propositions 7.4 and 7.6.

Choices of \([f_{ij}]\) in the numerical verification. Colorful cap: the entire range of admissible \([f_{ij}]\) subject to \(f_{12}^2+f_{13}^2+f_{23}^2=1\). Black and white dots: samples of \([f_{ij}]\) used in the numerical computation. Colors are only for visualization

Full size image

Choices of \(\{\omega _i\}\) in the numerical verification. triangular plate: the entire range of \(\{\omega _i\}\). dots: samples of \(\{\omega _i\}\). Portion enclosed by dashed line segments: \(\omega _1\leqslant \omega _2\leqslant \omega _3\) (Color figure online)

Full size image

Proposition 7.4

For \([f_{ij}]\) given by (5.6), the minimizers of (7.3) are (up to translation and reflection) the following with any \(l_k\geqslant 0\) satisfying \(\sum _{k=1}^nl_k=n\),

Proof

For (5.6), balls of type \(C\) do not engage in the interaction between charges. Similar to Proof of Proposition 7.1, we can prove that within any optimal arrangement, \(A\) and \(B\) must be alternating if \(C\) is ignored. Now let us exclude the following segment from optimal arrangements:

By translational invariance, we can assume that the above segment occupies the first 3 sites \(\{x_k\}_{k=1}^3\). Define the potential created by charges at all other sites as

Among the remaining balls centered at \(\{y_m\}_{m=4}^{3n}\), there are one more \(B\) than \(A\) , so \(V_2\) is piecewise quadratic in x with the coefficient of \(x^2\) being \(-\frac{1}{2}\), and thus is strictly concave on \([x_0,x_4]\). Therefore by swapping the \(A\) at \(x_2\) and the \(C\) at either \(x_1\) or \(x_3\), we can decrease U by the amount

Analogously, we can rule out other segments like \(C\)\(A\)\(B\)\(A\)\(C\) , \(C\)\(A\)\(B\)\(A\)\(B\)\(A\)\(C\) , etc. Lastly, let us verify that different choices of \(\{l_k\}\) yield the same U. Consider two \(A\)\(B\) dipoles separated by \(C\)

between which the potential energy is

and is independent of l. \(\square \)

Remark 7.5

In Proposition 7.4, if we penalize \(A\)\(C\) and \(B\)\(C\) interfaces with equal weights, then the minimizer is unique, i.e., \(A\)\(B\)\(\cdots \;\) \(A\)\(B\)\(C\)\(\cdots \;\) \(C\) . This is reminiscent of the results in Van Gennip and Peletier (2008, Theorems 4 and 5) where \(C\) forms only one macrodomain.

Proposition 7.6

For \([f_{ij}]\) given by (5.10) with nonpositive \(f_{12}\), \(f_{13}\) and \(f_{23}\), the arrangement \(A\)\(B\)\(C\)\(\cdots \)\(A\)\(B\)\(C\) minimizes (7.3).

Proof

By (5.10), \([f_{ij}]\) can be decomposed into three components, which are permutations of (5.6) and represent the interactions between sub-charges shown in Fig. 9. With the assumption \(f_{12},f_{13},f_{23}\leqslant 0\), all three components are simultaneously minimized by the cyclic arrangement \(A\)\(B\)\(C\)\(\cdots \)\(A\)\(B\)\(C\) , according to Proposition 7.4.   \(\square \)

C Numerical Computation of Free Energy in a 1-D Periodic Cell

We now offer some details of the numerical computation to seek for the 1-D global minimizers of the free energy (1.2). For each pattern, we obtain the optimal layer widths numerically, with the initial guess (for optimization) having uniform layer widths. We use , a constrained local minimization function in MATLAB, with the constraints being the volume constraints and nonnegativity of layer widths. For the input argument , we set , and to be \(10^{-6}\), which should be sufficient for our purposes. Since patterns are defined modulo translation (Ren and Wei 2003b, Definition 4.1), we can avoid some redundant computation. For further acceleration, we use MATLAB’s parallel tool to work on multiple (e.g., 24) patterns simultaneously.

We adopt a simple algorithm based on (4.3) to compute the long-range term of (1.2). Noticing that the Green’s function on [0, 1] with periodic boundary conditions is given by Ren and Wei (2003b, Equation (4.19))

we have \(\int _{y_1}^{y_2}\int _{x_1}^{x_2}G(x,y)\mathrm {d}{x}\mathrm {d}{y}=F(x_2\!-\!y_1)-F(x_1\!-\!y_1)-F(x_2\!-\!y_2)+F(x_1\!-\!y_2)\), where

Now, given a pattern and the positions of interfaces, the following algorithm returns the free energy J.

Although this algorithm is of complexity \(O\big (\) ⌃\(2\big )\), it is not a bottleneck compared to the exhaustive search (among all the patterns) of complexity \(O\big (2\) \(\big )\), since the total complexity is the product of the above two and that of (for constrained optimization in ).

D Analytic calculation of free energy of 1-D periodic patterns

For \(\Omega =[0,1]\) with periodic boundary conditions, the long-range term of (1.1) can be rewritten as

where \(v_j(x)=\int _0^1G(x,y)\,\big (u_j(y)-\omega _j\big )\mathrm {d}{y}\) (so that \(-v_j''=u_j-\omega _j\)) and \(w_j=-v_j'\,\). Denoting \({\vec {w}}=[w_1,w_2,w_3]^\mathrm{T}\), we can rewrite the above right-hand side as

Ren and Wei (2003b, Section 4), studied a local minimizer which is \(A\)\(B\)\(C\) identically repeating for n times. In that case, \({\vec {w}}\) is periodic with period 1/n, so one only needs to solve the following equation in one period in order to obtain the free energy

where \({\vec {\omega }}\) denotes \([\omega _1,\omega _2,\omega _3]^\mathrm{T}\), and \(\{\vec e_1,\vec e_2,\vec e_3\}\) forms the standard basis. In this way, we can obtain (3.1).

Analogously, for \(A\)\(B\)\(A\)\(C\) identically repeating for n times (with all the \(A\) layers having the same width), we can obtain (3.2) by solving the following equation

E Alternative Derivation of the Admissibility Conditions

As mentioned in Remark -(i), there is an alternative derivation of the conditions in Theorem 5.2 from the following three requirements:

• For \(\vec u={\vec {\omega }}\), the long-range term (4.2) attains zero;

• For any \(\vec u\) satisfying the incompressibility condition \(\vec u^\mathrm{T}{\vec {1}}=1\), the long-range term (4.2) is nonnegative;

• \([\gamma _{ij}]\) is symmetric.

Under Neumann or periodic boundary conditions, we have \(\int _{\Omega }G(\vec x,\vec y)\mathrm {d}{\vec x}=0\) for any \(\vec y\in \Omega \), so the first requirement is automatically satisfied and therefore does not lead to the condition \({\vec {\omega }}^\mathrm{T}[\gamma _{ij}]{\vec {\omega }}=0\). Moreover, the second requirement does not lead to the condition \([\gamma _{ij}]\succcurlyeq 0\) because of the incompressibility condition. However, as explained in Proposition 7.7, there are many equivalent choices of \([\gamma _{ij}]\), and one of them satisfies \([\gamma _{ij}]{\vec {\omega }}={\vec {0}}\) and \([\gamma _{ij}]\succcurlyeq 0\) as desired.

Proposition 7.7

Under the incompressibility condition, among all the \([\gamma _{ij}]\) fulfilling the above three requirements and yielding the same long-range term (4.2), there is a unique one satisfying \([\gamma _{ij}]\,{\vec {\omega }}={\vec {0}}\). Meanwhile, it also satisfies \([\gamma _{ij}]\succcurlyeq 0\).

Proof

Under the asssumption \(u_1+u_2+u_3=1\), we have \(\vec u=[u_1,u_2,u_3]^\mathrm{T}={\mathcal {A}}[u_1,u_2]^\mathrm{T}+[0,0,1]^\mathrm{T}\), where \({\mathcal {A}}\) is given by (7.5), and the second summand is constant. Since \(\int _{\Omega }G(\vec x,\vec y)\mathrm {d}{\vec x}=0\) for any \(\vec y\in \Omega \), we can rewrite (4.2) as

where \([{\tilde{\gamma }}_{ij}]={\mathcal {A}}^\mathrm{T}[\gamma _{ij}]{\mathcal {A}}\). To ensure that the above integral is nonnegative, we need to impose the condition \([{\tilde{\gamma }}_{ij}]\succcurlyeq 0\). (In fact, we can diagonalize \([{\tilde{\gamma }}_{ij}]\) into \(Q^\mathrm{T}\mathrm{diag}(\lambda _1,\lambda _2)\,Q\), and rewrite the above integral as a quadratic form like (4.4), from which it would be clear that \(\lambda _1\) and \(\lambda _2\) should be both nonnegative.)

By Lemma 7.8, there is a class of equivalent choices of \([\gamma _{ij}]\), but only one of them satisfies \([\gamma _{ij}]{\vec {\omega }}={\vec {0}}\). Such \([\gamma _{ij}]\) is given by (5.7) and is positive semi-definite since we have \([{\tilde{\gamma }}_{ij}]\succcurlyeq 0\). \(\square \)

Lemma 7.8

Define \(T:S_3\rightarrow S_2\) to be \(T(H)={\mathcal {A}}^\mathrm{T}H{\mathcal {A}}\), where \(S_m\) is the set of \(m\times m\) symmetric real matrices, and

then T is surjective. The kernel of T is \(\big \{{\vec {1}}\vec p^\mathrm{T}+\vec p{\vec {1}}^\mathrm{T}\;\big |\;\vec p\in {\mathbb {R}}^3\big \}\). Given any \({\vec {w}}=[w_1,w_2,w_3]^\mathrm{T}\in {\mathbb {R}}^3\) with \({\vec {w}}^\mathrm{T}{\vec {1}}=1\), in the quotient space \(S_3/\mathrm{ker}(T)\), the equivalence class of any \(H\in S_3\) has a unique representative \({\tilde{H}}\) satisfying \({\tilde{H}}\,{\vec {w}}={\vec {0}}\). This representative is given by \({\tilde{H}}={\mathcal {B}}^\mathrm{T}T(H){\mathcal {B}}\), where

Proof

We can take \(H=\begin{bmatrix}H_1 &{} 0\\ 0 &{} 0\end{bmatrix}\), where \(H_1\in S_2\), then \(T(H)=H_1\), so T is surjective. Since \({\vec {1}}^\mathrm{T}{\mathcal {A}}={\vec {0}}^\mathrm{T}\), based on the rank–nullity theorem we know \(\mathrm{ker}(T)=\big \{{\vec {1}}\vec p^\mathrm{T}+\vec p{\vec {1}}^\mathrm{T}\;\big |\;\vec p\in {\mathbb {R}}^3\big \}\). For any \(H\in S_3\), we have

where I is an identity matrix. According to the Sherman–Morrison formula, there is a unique \(\vec p\) such that the above right-hand side vanishes. We can verify that \({\tilde{H}}={\mathcal {B}}^\mathrm{T}T(H){\mathcal {B}}\) is the corresponding representative within the equivalence class of H, by checking \({\tilde{H}}{\vec {w}}={\vec {0}}\) and \(T({\tilde{H}})=T(H)\), which are clear from \({\mathcal {B}}{\vec {w}}={\vec {0}}\) and \({\mathcal {B}}{\mathcal {A}}=I\), respectively. \(\square \)

Remark 7.9

The results in this section can be generalized from \(\vec u\in {\mathbb {R}}^3\) to any dimension \({\mathbb {R}}^m\), by changing (7.5) into

and changing (7.6) into

where \(I_m\) is the \(m\times m\) identity matrix, and \({\vec {1}}\) is a vector whose components are all 1.

Reprints and permissions