[math(N)]개의 전자로 이루어진 분자의 파동함수를 표현하는데는 [math(4N)]개의 좌표가 필요하다.[1] 하지만 [math(N)]이 커질수록 파동함수는 우리가 필요한 정보 이상을 가지기도 하고 직관적으로 물리적 의미를 파악하기도 어려워진다. 그렇기에 파동함수 보다 더 적은 변수를 가지는 함수를 이용해 에너지나 다른 변수들의 기댓값을 구하려는 시도가 계속되었고, 그중 가장 유명한 이론이 전자 밀도를 이용한 밀도범함수 이론(Density Functional Theory, DFT)이다. 후술하겠지만 바닥 상태의 전자 밀도를 정확하게 안다면 신기하게도 계의 모양[2]이 하나로 정해지며 원론적으론 파동함수를 계산하지 않고도 모든 물리량의 기댓값을 계산할 수 있다! 물론 원론적으로만 그렇고 실제 전자 밀도를 실험 없이 바로 얻어낼 방법은 없으므로 이를 근사적으로 얻을 방법을 찾는 것이 DFT의 핵심이다.
어떻게 파동함수와 전자 확률 밀도를 연관지을 수 있을까? 우리가 찾고자 하는 양은 점 [math(\textbf{r})]에서 미소부피 [math(d\textbf{r}(=dxdydz))] 안에서 전자를 발견할 확률이다. 그리고 1번 전자가 [math(m_{s1})] 스핀을 가지고 좌표 [math(\textbf{r}_1)]에서 미소부피 [math(d\textbf{r}_1)]안에서 발견될 확률, 2번 전자가 [math(m_{s2})] 스핀을 가지고 좌표 [math(\textbf{r}_2)]에서 미소부피 [math(d\textbf{r}_2)]안에서 발견될 확률, [math(\cdots)], n번 전자가 [math(m_{sn})] 스핀을 가지고 좌표 [math(\textbf{r}_n)]에서 미소부피 [math(d\textbf{r}_n)]안에서 발견될 확률이
제1정리에선 바닥상태에서 축퇴가 없는 분자의 경우 바닥상태 전자 밀도 [math(\rho_0(\textbf{r}))][4]를 안다면 파동 함수와 모든 물리적 특징들이 하나로 결정됨을 보인다.[5] 이는 "[math(\rho_0(\textbf{r}))]이 해밀토니언을 결정한다."를 보이면 끝나는데 해밀토니언이 결정되면 파동함수와 에너지가 결정되고 파동함수를 통해 물리량 또한 결정되기 때문이다. 먼저 분자의 해밀토니언을 살펴보자. 보른-오펜하이머 근사를 적용시키면 전자의 해밀토니언은
와 같이 얻어진다. 위 식을 보면 [math(v(\textbf{r}))]과 전자 수 [math(n)]이 정해지면 해밀토니언은 자동으로 얻어짐을 알 수 있고 [math(\rho_0(\textbf{r}))]가 앞의 2개를 결정함을 보이면 된다. [math(\rho_0(\textbf{r}))]가 전자수를 결정한다는 것은 [math(n = \int \rho_0(\textbf{r}) d \textbf{r} )]을 통해 바로 알 수 있고 외부 퍼텐셜의 경우도 귀류법을 쉽게 통해 증명 가능하다. [math(\rho_0(\textbf{r}))]이 다른 외부 퍼텐셜[7] [math(v_a(\textbf{r}_i), v_b(\textbf{r}_i))]를 준다고 하면 해밀토니언도 정해지며 이를 [math(\mathcal{H}_a, \mathcal{H}_b)]라 쓰자. 바닥상태의 축퇴가 없다고 가정했으므로 각 해밀토니언은 단 1개의 정규화된 바닥상태 파동함수 [math(\psi_{0,a}, \psi_{0,b})]를 가질 것이다. 여기서 변분 원리를 쓰자.
증명은 간단한데 제1정리에 의해 [math(\rho_{\textup{trial}})]이 결정되면 [math(v_{\textup{trial}})]이 결정되고 [math(v_{\textup{trial}})]이 결정되면 [math(\rho_{\textup{trial}})]에 해당하는 [math(\psi_{\textup{trial}})]가 결정된다.[9] 그러면 변분 원리에 의해
호헨베르크-콘 이론을 통해 [math(E_0)]를 조건 [math(\int \rho \, d\textbf{r} - N = 0)] 하에서 [math(\rho)]에 대해 최소화 한다면 [math(\rho)]에 대한 방정식을 얻을 것이고 이를 풀면 이론적으론 [math(\rho)]를 파동함수 없이 구할 수 있다. 하지만 [math(\overline{T}[\rho])]와 [math(\overline{V}_{ee}[\rho])]의 꼴을 모르기에[11] 다른 방법을 고안해야 한다. 이 과정에서 핵심 역할을 하는 것이 콘-샴 방법(Kohn-Sham Method)이며 이론상으로 이 방법을 이용하면 정확한 결과값을 얻어낼 수 있다.[12]
[math(n)]개의 상호작용하지 않는 전자로 이루어져 있으며 모든 전자가 [math(v_s(\textbf{r}_i))]라는 외부 퍼텐셜 하에서 운동하고 있는 가상의 계를 생각해보자(reference system이라고 한다). 여기서 [math(v_s(\textbf{r}_i))]는 reference system의 바닥상태 전자 밀도 [math(\rho_s(\textbf{r}))]가 우리가 관심있어하는 계의 정확한 바닥상태 확률 밀도 [math(\rho_0(\textbf{r}))]와 동일하게 만들어주는 외부 퍼텐셜이다. 호헨베르크-콘 이론에서 증명했듯이 [math(\rho_s(\textbf{r}) = \rho_0(\textbf{r}))]로 전자 밀도를 설정하면 [math(v_s(\textbf{r}_i))]는 하나로 결정된다. 그럼 reference system의 해밀토니언은
가 된다. 여기서 문제는 [math(T_s[\rho])]의 꼴을 여전히 모른다는 것이다. 하지만 오비탈을 이용하면 표현할 수 있다. 이를 이용해 [math(\int \rho \, d\textbf{r} - N = 0)] 대신 [math(\int \theta^{KS*}_i\theta^{KS}_j \, d\textbf{r} - \delta_{ij} = 0)]라는 조건 하에서 [math(E)]를 최소화하는 문제로 바뀐다. 이는 라그랑주 승수법을 이용하면 다음과 같은 식[15]
[math(\displaystyle \begin{aligned} E_v = \sum\limits_{i} \langle i | T | i \rangle + \sum\limits_{i} \langle i | v | i \rangle + \frac{1}{2} \sum\limits_{i, j} \langle ij | ij \rangle + E_{xc} \end{aligned})]
가 된다. 이럼 [math(\Omega)]는
[math(\displaystyle \begin{aligned} \Omega = \sum\limits_{i} \langle i | T | i \rangle + \sum\limits_{i} \langle i | v | i \rangle + \frac{1}{2} \sum\limits_{i, j} \langle ij | ij \rangle + E_{xc} - \sum\limits_{i, j} \varepsilon_{ij} ( \langle i | j \rangle - \delta_{ij} ) \end{aligned})]
로 표현이 된다. 변분 계산 후 식 정리를 조금만 하면
[math(\displaystyle \begin{aligned} \delta \Omega &= \sum\limits_{i} ( \langle \delta i | T | i \rangle + \langle i | T | \delta i \rangle) + \sum\limits_{i} ( \langle \delta i | v | i \rangle + \langle i | v | \delta i \rangle ) \\ &+ \sum\limits_{i, j} ( \langle \delta ij | ij \rangle + \langle ij | \delta ij \rangle ) + \delta E_{xc} - \sum\limits_{i, j} ( \varepsilon_{ij} \langle \delta i | j \rangle - \varepsilon_{ji} \langle j | \delta i \rangle ) = 0 \end{aligned})]
콘-샴 방정식의 가장 큰 문제점은 [math(E_{xc})]와 [math(v_{xc})]의 모양을 모른다는 것이다. DFT 정확도의 핵심은 이 퍼텐셜을 실험치와 잘 부합하도록 근사하는 것이다. 그렇기에 정말 수많은 종류의 [math(E_{xc})]와 [math(v_{xc})]가 존재하며[17] 지금까지도 연구되고 있다. 아래는 그중 자주 쓰이는 일부만을 추려낸 것이다.
open-shell molecule들이나 해리 상태에 가까운 molecule들의 경우 LSDA(Local Spin Density Approximation)을 사용하는 것이 LDA보다 더 좋은 결과를 준다. HF의 UHF처럼 같은 콘-샴 오비탈에 들어가더라도 다른 모양을 줌으로써 오비탈의 유연성을 늘린다. [math(\alpha)] 스핀을 가지는 전자들의 밀도를 [math(\rho^{\alpha}(\textbf{r}))], [math(\beta)] 스핀을 가지는 전자들의 밀도를 [math(\rho^{\beta}(\textbf{r}))]라 하면 [math(E_{xc})]는 [math(\rho^{\alpha}(\textbf{r}))]와 [math(\rho^{\beta}(\textbf{r}))] 두개의 함수의 범함수가 된다.
이제 해야할 일은 LDA와 비슷하게 [math(\varepsilon_x(\rho^{\alpha}, \rho^{\beta}))]와 [math(\varepsilon_c(\rho^{\alpha}, \rho^{\beta}))]를 찾는 것이다. 보통 [math(\rho^{\alpha})], [math(\rho^{\beta})] 보다는 총 전자 밀도
LDA와 LSDA는 균일한 전자 가스로 부터 유도 되었기에 [math(\rho)]가 천천히 변할 때만 적절하다. 그렇기에 GGA(Generalized-Gradient Approximation)(Gradient-corrected functional이라고도 불린다.)에선 [math(\nabla \rho)] 항을 추가하여[21] 정확도를 높인다. 이를 식으로 표현하면
3\pi^2 \rho \right)^{1/3}\end{aligned} )]|| 여기서 쓰이는 매개변수 값은 두 방법 모두[math(\mu = 0.21951)]이고 PBE에선 [math(\kappa = 0.804)], revPBE에선 [math(\kappa = 1.245)]를 사용한다.
DFT에도 문제점들이 있는데 그중 하나가 장거리 전자 상관(Long-Range Electronic Correlation)을 고려하지 못한다는 것이다. 단거리 전자 상관(Short-Range Electronic Correlation)을 표현하는 것에만 신경을 쓰므로 공유 결합과 같은 것들은 잘 설명할지 몰라도 장거리 전자 상관이 중요한 전하 이동(Charge Transfer)나 분산력(Dispersion Force) 등을 설명하지 못한다는 단점이 있다. 이를 해결하기 위한 방법 중 하나론 RSH(Range-Separated Hybrid)이라 불리는 방법으로 다음과 같이 단거리에선 DFT Exchange를, 장거리에선 HF Exchange를 사용하는 것이다. 대표적인 예시인 CAM-B3LYP[48]에선
[1] 공간에 관한 좌표[math(3N)]개, 스핀에 관한 좌표 [math(N)]개[2] 핵의 위치, 핵의 전하량[3] 이는 [math(|\psi|^2)]가 두 입자의 변환에 대해 불변임을 이용하면 쉽게 증명 가능하다.[4] 0은 바닥상태임을 나타낸다.[5] 즉 앞서 언급한 것들은 [math(\rho(\textbf{r}))]의 범함수이다.[6] 이 양을 DFT에선 외부 퍼텐셜(external potential)이라 한다.[7] 상수 차이만 나는 퍼텐셜은 같은 퍼텐셜로 본다.[8] [math(v)]는 [math(E_0)]가 [math(v(\textbf{r}))]에 의존한다는 것을 강조하기 위해 붙었다.[9] 사실 이것은 [math(\rho_{\textup{trial}})]에 해당하는 반대칭 파동함수를 주는 [math(v_{\textup{trial}})]가 존재할때만 성립한다. 이러한 [math(\rho_{\textup{trial}})]를 [math(v)]-representable하다고 하며 모든 [math(\rho_{\textup{trial}})]가 [math(v)]-representable하지는 않다. 다만 크게 걱정할 필요는 없는게 [math(v)]-representable 조건을 약화시킬 수 있는데 이는 Density-Finctional Theory of Atoms and Molecules를 참고하여라.[10] 이 또한 Density-Finctional Theory of Atoms and Molecules를 참고하여라.[11] 물론 이 항을 근사적으로 나타낸다면 이 방법 또한 [math(\rho)]를 구하는데 쓸 수 있다. 이를 Orbital-free density functional theory라 한다. 하지만 정확도 문제로 자주 쓰이는 방법은 아니다.[12] 물론 이론상으로 그렇다는 것이고 실제로는 근사가 필요하다.[13] 이제부터 첨자 0은 빼도록 하겠다.[14]전기 퍼텐셜 참조[15] 아래에서 나오는 표기법들은 하트리-포크 방법 참고[16] 잘 쓰이는 방법은 아니다.[17]참고[18]오일러-라그랑주 방정식 참조[19] P. A. M. Dirac (1930) P. Camb. Philos. Soc. 26, pp. 376.[20] S. H. Vosko, L. Wilk, and M. Nusair (1980) Can. J. Phys. 58, pp. 1200.[21] [math(\nabla^2 \rho)]까지 추가하기도 한다.[22] A. D. Becke (1988) Phys. Rev. A 38, pp. 3098.[23] J. P. Perdew, K. Burke, and M. Ernzerhof (1996) Phys. Rev. Lett. 77, pp. 3865.[24] Y. Zhang and W. Yang (1998) Phys. Rev. Lett. 80, pp. 890.[25] J. P. Perdew, K. Burke, and M. Ernzerhof (1996) Phys. Rev. Lett. 77, pp. 3865.[26] C. Lee, W. Yang, and R. G. Parr (1988) Phys. Rev. B 37, pp. 785.[27] J. Tao, J. P. Perdew, V. N. Staroverov, and G. E. Scuseria (2003) Phys. Rev. Lett. 91, pp. 146401.[28] J. P. Perdew, A. Ruzsinszky, G. I. Csonka, L. A. Constantin, and J. Sun (2009) Phys. Rev. Lett. 103, pp. 026403.[29] Y. Zhao and D. G. Truhlar (2006) J. Chem. Phys. 125, pp. 194101.[30]전이 금속, 무기화합물, 유기금속화합물 계산에 주로 쓰인다.[31] A. D. Becke (1993) J. Chem. Phys. 98, pp. 5648.[32] P. J. Stephens, F. J. Devlin, C. F. Chabolowski, and M. J. Frisch (1994) J. Phys. Chem. 98, pp. 11623.[33] C. Adamo and V. Barone (1999) J. Chem. Phys. 110, pp. 6158.[34] A. D. Becke (1997) J. Chem. Phys. 107, pp. 8554.[35] Y. Zhao and D. G. Truhlar (2008) Theor. Chem. Acc. 120, pp. 215.[36] 주족 원소 열화학, 동역학, 비공유결합성 상호작용 계산에선 매우 좋은 성능을 지니지만 전이 금속 열화학이나 유기금속화합물 계산엔 사용할 수 없다.[37] Y. Zhao and D. G. Truhlar (2007) J. Chem. Theory Comput. 4, pp. 1849.[38] 주족 원소 열화학, 동역학, 비공유결합성 상호작용 계산에 쓰인다.[39] V. N. Staroverov, G. E. Scuseria, J. Tao, and J. P. Perdew (2003) J. Chem. Phys. 119, pp. 12129.[40] G. I. Csonka, J. P. Perdew, and A. Ruzsinszky (2010) J. Chem. Theory Comput. 6, pp. 3688.[41] 한 예시로 MP2의 이차 보정에너지를 추가한다.[42] S. Grimme (2006) J. Chem. Phys. 124, pp. 034108.[43] Y. Zhang, X. Xu, and W. A. Goddard III (2009) Proc. Natl. Acad. Sci. USA 106, pp. 4963.[44] 특이하게 [math(E^{\textup{MP2}}_{c})]의 계산엔 B3LYP의 오비탈을 사용한다.[45] I. Y. Zhang, X. Xin, Y. Jung, and W. A. Goddard III (2011) Proc. Natl. Acad. Sci. USA 108, pp. 19896.[46] S. Kozuch and J. M. L. Martin (2013) J. Comput. Chem. 34, pp. 2327.[47] L. Goerigk and S. Grimme (2011) J. Chem. Theory Comput. 7, pp. 291.[48] T. Yanai, D. P. Tew, and N. C. Handy (2004) Chem. Phys. Lett. 393, pp. 51.[49] N. Mardirossian and M. Head-Gordon (2014) Phys. Chem. Chem. Phys. 16, pp. 9904.[50] Y.-S. Lin, G.-D. Li, S.-P. Mao, and J.-D. Chai (2013) J. Chem. Theory Comput. 9, pp. 263.[51] J.-D. Chai and M. Head-Gordon (2008b) Phys. Chem. Chem. Phys. 10, pp. 6615.[52] N. Mardirossian and M. Head-Gordon (2016) J. Chem. Phys. 144, pp. 214110.[53] 이론상 그렇다는 것이고 실제론 [math(v_{xc})]를 정확하게 찾는 것은 불가능하다.