混合の自由エネルギー変化をプロットしてみたよ

# 混ざるか混ざらないか
液体同士を混ぜ合わせた時に混ざるものと混ざらないものがある。透き通った紅茶にミルクを入れると全体が均一に濁って美味しいミルクティーになる。かたや、水に油を落としてよく振ってみても二層に分かれる。両方とも液体を混ぜ合わせるのは同じなのに、簡単に混ざるものと全く混ざらないものがある。この違いはなんなのだろうか。

# 混合の自由エネルギー変化

熱力学的にいうと、定圧状態において反応が自発的に進行するかどうかは、系の始状態と終状態のギブスの自由エネルギーを比較して、減少していれば自発的に進行すると言える。今回の場合は液体１、液体２を混ぜる前と混ぜた後の自由エネルギーの変化を見れば良い。

自由エネルギー$ G$は次のように定義される。

$ G = H - TS$

ここで、混合によるエントロピー変化$\Delta S_{mix}$とエンタルピー変化$\Delta H_{mix}$は、次のようになる。導出は省略します。

$ \Delta S_{mix} = -R [n_1 \log(x_1) + n_2 \log(x_2)]$

$ \Delta H_{mix} = R T \chi \frac{n_1 n_2}{n_1 + n_2}$

$n_1$, $x_1$は液体1の物質量とモル分率、$n_2$, $x_2$は液体2の物質量とモル分率であり、$ \chi$は溶解パラメーターを表す。溶解パラメーターは無次元の値である。この二つをギブスの自由エネルギーの式に代入すると混合の自由エネルギー変化$ \Delta G_{mix}$が得られる。

$ \Delta G_{mix} = R T [n_1 \log(x_1) + n_2 \log(x_2) + \chi \frac{n_1 n_2}{n_1 + n_2} ]$



上記の理論は高分子化学においてFlory–Huggins理論と呼ばれ非常によく用いられている。今回は低分子同士の混合について考えているのでFlory–Huggins理論における高分子鎖セグメント長が1の場合についての式を示してる。導出等の詳細はFlory–Huggins理論と検索していただき、そちらをご参照ください。

# プロットしよう
今回は混合について少し議論を深めるために三種類のプロットを行う。

## 混合の自由エネルギー変化の概形の確認
まずは定数を固定して自由エネルギー変化の概形を確認しよう。

```Python
import matplotlib.pyplot as plt
import numpy as np

R = 8.314 # [J/K mol]
T = 300 # [K]

kai = 3

n1 = np.linspace(0,1,10000)
n2 = 1 - n1

x1 = n1/(n1+n2)
x2 = n2/(n1+n2)

fig = plt.figure()
ax = fig.add_subplot(111)

f0 = 0*x1
plt.plot(x1,f0,'--', color='black')

F = R*T*(n1*np.log(x1) + n2*np.log(x2) + kai*n1*n2/(n1 + n2))
F1 = R*T*n1*np.log(x1)
F2 = R*T*n2*np.log(x2)
F3 = R*T*kai*n1*n2/(n1 + n2)

plt.plot(x1,F)
plt.plot(x1,F1,'--', color = 'red')
plt.plot(x1,F2,'--', color = 'orange')
plt.plot(x1,F3,'--', color = 'green')

ax.set_xlabel('x$_1$', fontsize=20, color='black')
ax.set_ylabel('free energy', fontsize=20, color='black')
```
今回、液体1を0 molから1 molまで変化させている。液体2の物質量は$n_2 = 1 - n_1$である。温度を300 K、溶解パラメーターを3に設定してみた。
実行結果を以下に示そう。

![mix1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2804113/2e6e36f3-d036-5680-daa5-7ab9032cb8eb.png)

いい感じだと思います。
縦軸が自由エネルギーで横軸が液体1の分率である。青実線が$ \Delta G_{mix}$、赤点線が$ R T n_1 \log(x_1)$、黄色点線が$ R T n_2 \log(x_2)$、最後に緑点線が$ \Delta H_{mix}$である。ここから、$ \Delta H_{mix}$は上に凸の関数、$ \Delta S_{mix}$は下に凸の関数でそれらの足し合わせが$ \Delta G_{mix}$となる。$ \Delta G_{mix}$が負の時二つの溶媒は自発的に混ざり合うので、今回の条件では$ x_1$がおよそ0から0.2、0.8から1の時にのみ混ざり合うことがわかる。

## 溶解パラメーターを変化させてみる
今度は温度を一定にして溶解パラメーターを変化させてみた図をプロットしてみよう。作図のプログラムは、上のプログラムとほとんど一緒なので重なる部分も多いですが全部書いておきます。

```Python
import matplotlib.pyplot as plt
import numpy as np

R = 8.314 #[J/K mol]
T = 300 #[K]

kai = np.linspace(-1,5,10)

n1 = np.linspace(0,1,10000)
n2 = 1 - n1

x1 = n1/(n1+n2)
x2 = n2/(n1+n2)

fig = plt.figure()
ax = fig.add_subplot(111)

f0 = 0*x1
plt.plot(x1,f0,'--', color='black')

for i in range(10):
    F = R*T*(n1*np.log(x1) + n2*np.log(x2) + kai[i]*n1*n2/(n1 + n2))
    
    plt.plot(x1,F)

ax.set_xlabel('x$_1$', fontsize=20, color='black')
ax.set_ylabel('free energy', fontsize=20, color='black')
```
今回は$ \chi$を-1から5まで振ってみたがどうだろうか。実行結果を以下に示す。

![mix2.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2804113/169ed19e-808f-33e9-12d8-2300514416da.png)

いい感じにできましたね。
溶解パラメーターが小さいときは$ \Delta G_{mix}$が全範囲で負の値を示し、$ \chi$が大きくなるにつれて、中心が正の値となり、最後にはほぼ全範囲で正の値となる。$ \chi$が小さいときはどんな分率においても混ざり合い、大きくなると$ x_1 = 5$付近から混ざらない状態が出てきて、溶解パラメーターが大きいときはほとんど全ての分率で混ざり合わないことに対応する。溶解パラメーターの物理的意味は液体1と液体2の親和性みたいなものを示しており、$ \chi$が小さいと親和性が高いといえる。

## 温度を変化させてみる
今度は$ \chi = 3$に固定して温度を振ってみよう。この条件のイメージとしては、同じ溶液同士を温度を変えて混ぜてみるみたいな感じです。温度変えると大抵の場合溶解度って大きくなると思うので$ \Delta G_{mix} < 0$の範囲が増えるのではないでしょうか。再度ほとんど同じプログラムなので、冗長ではありますが全文載せておきます。

```Python
import matplotlib.pyplot as plt
import numpy as np

R = 8.314 # [J/K mol]
T =  np.linspace(50,500,10)#[K]

kai = 3

n1 = np.linspace(0,1,10000)
n2 = 1 - n1

x1 = n1/(n1+n2)
x2 = n2/(n1+n2)

fig = plt.figure()
ax = fig.add_subplot(111)

f0 = 0*x1
plt.plot(x1,f0,'--', color='black')

for i in range(10):
    F = R*T[i]*(n1*np.log(x1) + n2*np.log(x2) + kai*n1*n2/(n1 + n2))
    
    plt.plot(x1,F)

ax.set_xlabel('x$_1$', fontsize=20, color='black')
ax.set_ylabel('free energy', fontsize=20, color='black')
```

今度は$ \chi$を3に固定し温度を50 Kから500 Kまで変化させてみた。次に実行結果を載せます。

![mix3.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2804113/72e3f05b-60c4-b6a1-374d-5dbe84b8b1fb.png)


ありゃりゃ、ちょっと意外だけど面白いのでいいのではないでしょうか。
温度変えても溶ける分率は変わらないという結果になりました。もちろんモデル上計算値がこうなるというだけで、現実世界でも温度を変えても溶ける分率が変わらないということではないと思います。実際、温度を上げると大抵の物質の溶解度は上昇するので違うのではないでしょうか。直感と計算結果が異なるのでモデルの限界かもしれないと疑ってみてどんな原因があるか考えてみます。

まず、Flory–Huggins理論は格子モデルに基づいているので溶解時の体積変化が考慮できない理論となっております。温度変化による体積変化もまた影響があるかもしれません。また、現実の世界では物質同士の親和性も温度に依存して変化しそうなので、$\chi$を固定して$T$だけ変化させるっていうのも現実と乖離しているのかもしれません。

直感と違うというだけで、実際測定してみると理論通りになっているのかもしれません。時間がある時に、調べてみようと思います。でも、実際モデルが違うなら、モデルの限界をプロットから検討できたということになって面白いなあと思います。やっぱ図示してみるのって大切なんだなあと思いました。ありがとうございました。