[Physics] Flory-Huggins ternary phase diagram with a neutral component

phase diagramphase-transitionpolymersstatistical mechanicsthermodynamics

I am searching the literature for the Flory-Huggins phase diagram with the following components : polymer, solvent, and a third component that does not interact with the other components (just entropy effects). It must have been done but I can't find the special case in which the third component is not interacting. I am interested in the entropy effect of the third component.

Thanks for any help.

Best,
Jd

Best Answer

Introduction

In this answer, I first tackle the simpler problem of a mixture with three types of molecules identical in size. This has relevance to the original question in that the third component does not directly interact with the others, yet its introduction changes the phase diagram through entropic effects.

Small molecules

Take the Flory-Huggins free energy of the form $$\frac{\mathcal{F}}{k_BT} = \chi\varphi_\text{A}\varphi_\text{B} + \varphi_\text{A}\ln\varphi_\text{A} + \varphi_\text{B}\ln\varphi_\text{B} + \varphi_\text{C}\ln\varphi_\text{C}$$

Here the molecule types A, B have an enthalpic interaction between them characterized by the Flory-Huggins parameter $\chi$, and the component C only appears in the entropy of mixing -term. On that note, the present free energy can be obtained by a three state (Potts) Ising model via the Bragg-Williams mean field approximation. The discussion here is qualitative, but for more quantitative theories one ought to use semiempirical equations of state, or if the equations come from toy models, more advanced approximative methods. (Bethe) Guggenheim's quasichemical approximation comes to mind as perhaps the simplest of such methods, and can be numerically evaluated with relatively little extra effort.

Taking the opposite point of view, one might want to try and make simplifying approximations to the free energy, for even for the present simple case, its minima cannot be evaluated without the use of numerical methods. The (Landau) $\varphi^4$ theory, where you Taylor expand the free energy w.r.t. the order parameter to the fourth order is perhaps the most popular approximation. The odd powers drop out in the expansion due to the symmetry of the system, and thus $\varphi^4$ is the first polynomial degree where you actually can have two phases.

To actually go about solving the case at hand. one could assume $\varphi_\text{C}$ constant and then write $$\frac{\mathcal{F}}{k_BT} = \chi\varphi_\text{A}(c-\varphi_\text{A}) + \varphi_\text{A}\ln\varphi_\text{A} + (c-\varphi_\text{A})\ln(c-\varphi_\text{A}) + \varphi_\text{C}\ln\varphi_\text{C}$$ where for notational convenience I wrote $c = 1-\varphi_\text{C}$. The only thing that is then truly different from the two-component case is that the critical value for the FH parameter is $\chi = 2/c$ rather than $2$. Now obviously in holding $\varphi_\text{C}$ constant, we have assumed that the tielines should be such that the C component is in the same proportion in the two phases. This is indeed the case due to the symmetry of the system.

I have graphed the ternary phase diagram it below for $\chi = 3$ (horizontal tielines). Phase diagram

General case

For polymers the Flory-Huggins free energy is slightly more complicated, $$\frac{\mathcal{F}}{k_BT} = \chi\varphi_\text{A}\varphi_\text{B} + \varphi_\text{A}\ln\varphi_\text{A} + \frac{\varphi_\text{B}}{N}\ln\varphi_\text{B} + \varphi_\text{C}\ln\varphi_\text{C}$$

where $N$ is the size of the polymer (often one makes approximations where $N\gg1$). This simple addition considerably changes how the system behaves, for the symmetry between A and B is then fundamentally broken.

The analysis is also much more involved (at least my numerical implementation is) than the Small molecules case for two reasons: The system is not perfectly symmetric so it does not phase separate to two different free energy minima of the uniformly mixed system. Rather, one has to find the points with the same chemical potential, and then find the one that produces the wanted composition with the lowest free energy. Second, the tielines need not be horizontal, so there is another dimension with respect to which to minimize.

I have provided my code in Appendix A, with which I have produced the graph below with $N = 2$, $\chi = 3$. For your question I suppose the most important aspect of this is that the tielines are indeed not horizontal. Phase diagram

Appendix A: General case

The following "Python" code hack does the job. It has some numerical issues when the arguments in the logarithms go near zero and such. This could probably be alleviated by Taylor expansions near the borders. As it stands, the code is uncommented and horribly ugly. If I have time, I'll clean it down and explain where some of the horrible polynomials etc. come from (handling the domain borders is a mess). Do not trust the answers this spits at you, I have yet to properly verify that it works, but I thought you might want to have a look at the code, however dysfunctional, sooner rather than later.

from scipy.optimize import *
# Size of polymer
N = 1.
# chi parameter
w = 3.

# function defs
def freeenfun(x, w):
  phia = x[0]; phib = x[1]; phic = 1.-phia-phib
  return w*phia*phib + phia*log(phia) + phib/N*log(phib) + phic*log(phic)
def freeender(x, w):
  phia = x[0]; phib = x[1]; phic = 1.-phia-phib
  fda  = w*phib + log(phia/phic)
  fdb  = w*phia + log(phib)/N + 1./N - log(phic) - 1.
  return fda, fdb
def freeenderdir(t, x, y, w, lvl):
  phia = x[0] + t*cos(y); phib = x[1] + t*sin(y); phic = 1.-phia-phib
  fda, fdb = freeender(array([phia, phib]), w)
  fddir = fda*cos(y) + fdb*sin(y)
  return fddir - lvl
def freeenfun2(t, x, y, w):
  phia = x[0] + t*cos(y); phib = x[1] + t*sin(y); phic = 1.-phia-phib
  return freeenfun(array([phia, phib]), w)
def frompatophi(t, x, y):
  phia = x[0] + t*cos(y); phib = x[1] + t*sin(y); phic = 1.-phia-phib
  return phia, phib, phic

def getonepair(y, x, w): 
  y = y[0]
  a = cos(y)
  b = sin(y)
  bpol = roots(array([ (2*a**2*b**3 + 2*a**3*b**2)*w*N,
  ((4*a**2*b**2 + 2*a**3*b)*w*x[1] + (2*a*b**3 + 4*a**2*b**2)*w*x[0] - 2*a**2*b**2*w - a*b**3 - a**2*b**2)*N + a*b**3 + a**2*b**2,
  (2*a**2*b*w*x[1]**2 + ((4*a*b**2 + 4*a**2*b)*w*x[0] - 2*a**2*b*w - a*b**2)*x[1] + 2*a*b**2*w*x[0]**2 + (-2*a*b**2*w - b**3 - 2*a*b**2)*x[0] - a**2*b)*N + a*b**2*x[1] + (b**3 + 2*a*b**2)*x[0] - a*b**2,
  ((2*a*b*w*x[0] + a**2)*x[1]**2 + (2*a*b*w*x[0]**2 + (-2*a*b*w - b**2 - 2*a*b)*x[0] - a**2)*x[1])*N + b**2*x[0]*x[1] + b**2*x[0]**2 - b**2*x[0] ]))[1:]
  if any(iscomplex(bpol)):
    return array([nan, nan])
  lowerb = min(bpol)
  upperb = max(bpol)
  fdmax = freeenderdir(lowerb, x, y, w, 0)
  fdmin = freeenderdir(upperb, x, y, w, 0)
  hh = x[0]*tan(pi-y)
  ww = x[1]*tan(y-pi/2)
  if hh + x[1] < 1.:
    tmax = sqrt(x[0]**2 + hh**2)
  else:
    tmax = sin(pi/4.)/sin(3.*pi/4.-y)*(1-x[1]-x[0])
  if ww + x[0] < 1.:
    tmin = -sqrt(x[1]**2 + ww**2)
  else:
    tmin = -sin(pi/4.)/sin(y-3.*pi/4.)*(1-x[0]-x[1])
  def getonepair2(fd):
    ma = brentq(freeenderdir, tmin+1e-12, lowerb, args=(x, y, w, fd))
    mb = brentq(freeenderdir, upperb, tmax-1e-12, args=(x, y, w, fd))
    return array([ma, mb])
  def minfd(fd):
    onepair = getonepair2(fd)
    ypai0 = freeenfun2(onepair[0], x, y, w)
    ypai1 = freeenfun2(onepair[1], x, y, w)
    numd = (ypai1-ypai0)/diff(onepair)
    return (freeenderdir(onepair[0], x, y, w, 0) - numd)**2
  minres = minimize(minfd, (fdmax+fdmin)/2.)#, method='SLSQP', bounds=((fdmin, fdmax),))
  return getonepair2(minres.x)

def minang(x, w): 
  def minfun(y):
    onepair = getonepair(y, x, w)
    a = -onepair[0]/(onepair[1]-onepair[0])
    ff = (1-a)*freeenfun2(onepair[0], x, y, w) + a*freeenfun2(onepair[1], x, y, w)
    return ff
  minres = minimize(minfun, 135./180*pi)
  return minres.x, getonepair(minres.x, x, w)

aa = linspace(0.2, 0.4999, 256)
phia = ones((2, len(aa)))*0
phib = zeros_like(phia)
phic = zeros_like(phia)
for i in range(len(aa)):
  x = array([1, 1])*aa[i]; 
  y, tt = minang(x, w)
  phia[:, i], phib[:, i], phic[:, i] = frompatophi(tt, x, y)

kk = vstack((phia[0, :], phib[0, :], phic[0, :], phia[1, :], phib[1, :], phic[1, :])).transpose()
kk = kk[-isnan(kk[:, 0]), :]
savetxt('tern.dat', kk, header='A_a A_b A_c B_a B_b C_c', comments='')

Appendix B: Plotting ternary axes

This is the code snippet I use for my graphs.

\documentclass[class=minimal,border=0pt]{standalone}
\usepackage{pgfplots}
\usepgfplotslibrary{ternary}
\begin{document}
\begin{tikzpicture}
\begin{ternaryaxis}[
    ylabel=A,
    zlabel=B,
    xlabel=C
]
    \addplot3[smooth,
        tieline={each nth tie=10},
        fill=blue!10,
    ]
    table [y=A_a,z=A_b,x=A_c] {tern.dat};
\end{ternaryaxis}
\end{tikzpicture}

\end{document}
Related Question