-
-
Save ohno/d63d79ca816a3614b201bc5728399fb7 to your computer and use it in GitHub Desktop.
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| { | |
| "cells": [ | |
| { | |
| "cell_type": "markdown", | |
| "id": "ee76f64e", | |
| "metadata": {}, | |
| "source": [ | |
| "# モンテカルロ積分を用いたヘリウム原子の摂動計算\n", | |
| "\n", | |
| "ここではヘリウム原子の電子状態計算を題材として摂動法を学びます. 解析的な計算を避けて気軽に計算できるよう波動関数には[Antique.jl](https://ohno.github.io/Antique.jl/dev/HydrogenAtom/)を使用し, 積分には[MonteCarloIntegration.jl](https://github.com/ranjanan/MonteCarloIntegration.jl)を使用します. なお, 摂動論のちゃんとした説明は[東京工業大学の武藤先生の講義ノート](http://www.th.phys.titech.ac.jp/~muto/lectures/QMII09/QMII09_chap14.pdf)を参照してください." | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "id": "e46f90de", | |
| "metadata": {}, | |
| "source": [ | |
| "## パッケージの使い方\n", | |
| "\n", | |
| "ヘリウム原子$\\mathrm{He}$を扱う前に, まずはヘリウム原子イオン$\\mathrm{He}^+$を題材としてパッケージの使い方を学びましょう. $\\mathrm{He}^+$は$+2$の電荷をもった核と$-1$の電荷をもった電子が1つずつ集まってできている, いわゆる水素様原子（すいそようげんし）です. 水素様原子の固有値・固有関数は[Antique.jl](https://ohno.github.io/Antique.jl/dev/)で利用できます. まず下記を実行してインストールしましょう." | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 1, | |
| "id": "512ff3a3", | |
| "metadata": {}, | |
| "outputs": [], | |
| "source": [ | |
| "# using Pkg; Pkg.add(\"Antique\")" | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "id": "68c8edd4", | |
| "metadata": {}, | |
| "source": [ | |
| "$\\mathrm{He}^+$は原子番号（核の電荷）が$Z=2$なので, [こちら](https://ohno.github.io/Antique.jl/dev/HydrogenAtom/#Usage-and-Examples)を例にして`Z=2`を指定します." | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 2, | |
| "id": "cbfa3a37", | |
| "metadata": {}, | |
| "outputs": [ | |
| { | |
| "data": { | |
| "text/plain": [ | |
| "Base.Meta.HydrogenAtom15725040630537206220" | |
| ] | |
| }, | |
| "execution_count": 2, | |
| "metadata": {}, | |
| "output_type": "execute_result" | |
| } | |
| ], | |
| "source": [ | |
| "using Antique\n", | |
| "He⁺ = antique(:HydrogenAtom, Z=2, Eₕ=1.0, a₀=1.0, mₑ=1.0, ℏ=1.0)" | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "id": "cc391f6e", | |
| "metadata": {}, | |
| "source": [ | |
| "核の質量が十分に大きいときの水素様原子では, 主量子数$n \\geq 1$に対する固有値は\n", | |
| "\n", | |
| "$$\n", | |
| "E_n = -\\frac{Z^2}{2 n^2} E_{\\mathrm{h}}\n", | |
| "$$\n", | |
| "\n", | |
| "で与えられるので, $\\mathrm{He}^+$の基底状態, つまり$Z=2, n=1$では\n", | |
| "\n", | |
| "$$\n", | |
| "E_1 = -\\frac{2^2}{2 1^2} E_{\\mathrm{h}} = -2 E_{\\mathrm{h}}\n", | |
| "$$\n", | |
| "\n", | |
| "となります. これは先ほど定義したモジュール`He⁺`の関数`E()`に`n=1`を指定して次のように計算できます." | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 3, | |
| "id": "4bb06883", | |
| "metadata": {}, | |
| "outputs": [ | |
| { | |
| "data": { | |
| "text/plain": [ | |
| "-2.0" | |
| ] | |
| }, | |
| "execution_count": 3, | |
| "metadata": {}, | |
| "output_type": "execute_result" | |
| } | |
| ], | |
| "source": [ | |
| "He⁺.E(n=1)" | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "id": "add29dc4", | |
| "metadata": {}, | |
| "source": [ | |
| "エネルギーの単位について, [2018年CODATA推奨値](https://physics.nist.gov/cgi-bin/cuu/Value?hrev)の換算係数を用いると次のような値に換算できます." | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 4, | |
| "id": "d6f79d91", | |
| "metadata": {}, | |
| "outputs": [ | |
| { | |
| "name": "stdout", | |
| "output_type": "stream", | |
| "text": [ | |
| "E = -2.0 Eₕ\n", | |
| " = -54.422772491976 eV\n" | |
| ] | |
| } | |
| ], | |
| "source": [ | |
| "Eₕ2eV = 27.211386245988 # eV # https://physics.nist.gov/cgi-bin/cuu/Value?hrev\n", | |
| "println(\"E = \", He⁺.E(n=1), \" Eₕ\")\n", | |
| "println(\" = \", He⁺.E(n=1) * Eₕ2eV, \" eV\")" | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "id": "56c3ffbf", | |
| "metadata": {}, | |
| "source": [ | |
| "次に, 波動関数の規格化条件を確認しましょう. 積分には[MonteCarloIntegration.jl](https://github.com/ranjanan/MonteCarloIntegration.jl)を使用するので次のコマンドを実行してインストールします." | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 5, | |
| "id": "9814b4bf", | |
| "metadata": {}, | |
| "outputs": [], | |
| "source": [ | |
| "# using Pkg; Pkg.add(\"MonteCarloIntegration\")" | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "id": "2e15be27", | |
| "metadata": {}, | |
| "source": [ | |
| "このパッケージで次の積分を計算する例です.\n", | |
| "\n", | |
| "$$\n", | |
| "\\int |\\psi_n(r,\\theta,\\varphi)|^2 r^2 \\sin\\theta~\\mathrm{d}r \\mathrm{d}\\theta \\mathrm{d}\\varphi = 1\n", | |
| "$$" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 6, | |
| "id": "881775ef", | |
| "metadata": {}, | |
| "outputs": [ | |
| { | |
| "name": "stdout", | |
| "output_type": "stream", | |
| "text": [ | |
| "<ψ|ψ> = 0.9983502045775035 ± 0.0007877199112562297\n" | |
| ] | |
| } | |
| ], | |
| "source": [ | |
| "using MonteCarloIntegration\n", | |
| "res = vegas(r -> abs(He⁺.ψ(r...))^2 * r[1]^2 * sin(r[2]), [0,0,0], [100,π,2π], atol=1e-5)\n", | |
| "println(\"<ψ|ψ> = \", res.integral_estimate, \" ± \", res.standard_deviation)" | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "id": "fe1d8043", | |
| "metadata": {}, | |
| "source": [ | |
| "おおよそ1が得られましたので, 正しく計算できると考えられます. プログラムについて補足します. まず`r[1]`は$r$, `r[2]`は$\\theta$, `r[3]`は$\\varphi$に対応します. これだけ説明すれば後は察せられると思いますが, 本来ならば$0\\leq r \\lt \\infty, 0\\leq \\theta \\lt \\pi, 0\\leq \\varphi \\lt 2\\pi$を計算するところを, $r$に関しては無限遠までサンプリングすることはできないので適当に$r=100$程度までで打ち切っています. また, `He⁺.ψ(r...)`は[splat](https://qiita.com/tenfu2tea/items/5c0ecb10dfdd2002e262#%E5%AE%9F%E5%BC%95%E6%95%B0%E3%81%AE%E3%83%AA%E3%82%B9%E3%83%88%E3%82%92%E8%A4%87%E6%95%B0%E3%81%AB%E5%B1%95%E9%96%8B%E3%81%99%E3%82%8B-splat)記法を用いた`He⁺.ψ(r[1],r[2],r[3])`の略記です." | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "id": "f12c0a50", | |
| "metadata": {}, | |
| "source": [ | |
| "## ハミルトニアン\n", | |
| "\n", | |
| "ヘリウム原子を2つの電子と1つの原子核（アルファ粒子＝中性子2つと陽子2つ）からできた3体系として扱います. 電子と核はいずれも質点と考え, 相互作用はクーロン力のみとすれば, 非相対論的なハミルトニアンは次式で書き表されます.\n", | |
| "\n", | |
| "$$\n", | |
| "\\hat{H}_\\mathrm{tot}\n", | |
| "=\n", | |
| "- \\frac{\\hbar^2}{2 m_\\mathrm{\\alpha}} \\nabla_{\\alpha}^2\n", | |
| "- \\frac{\\hbar^2}{2 m_\\mathrm{e}} \\nabla_{1}^2\n", | |
| "- \\frac{\\hbar^2}{2 m_\\mathrm{e}} \\nabla_{2}^2\n", | |
| "- \\frac{2 e^2}{4 \\pi \\varepsilon r_1}\n", | |
| "- \\frac{2 e^2}{4 \\pi \\varepsilon r_2}\n", | |
| "+ \\frac{ e^2}{4 \\pi \\varepsilon r_{12}}\n", | |
| "$$\n", | |
| "\n", | |
| "これを原子単位（詳細は[IUPAC Greenbook 3.2.9](https://unit.aist.go.jp/nmij/public/report/others/pdf/iupac_green_book_jp.pdf)を参照）で無次元化すると\n", | |
| "\n", | |
| "$$\n", | |
| "\\hat{H}_\\mathrm{tot}^\\ast =\n", | |
| "- \\frac{1}{2 m_\\mathrm{\\alpha}/m_\\mathrm{e}} \\nabla_{\\alpha}^{2\\ast}\n", | |
| "- \\frac{1}{2} \\nabla_{1}^{2\\ast}\n", | |
| "- \\frac{1}{2} \\nabla_{2}^{2\\ast}\n", | |
| "- \\frac{2}{r_1^\\ast}\n", | |
| "- \\frac{2}{r_2^\\ast}\n", | |
| "+ \\frac{1}{r_{12}^\\ast}\n", | |
| "$$\n", | |
| "\n", | |
| "となります. ここで, ハミルトニアンはハートリーエネルギー$E_\\mathrm{h}$を用いて無次元化されており, $\\hat{H}_\\mathrm{tot} / E_\\mathrm{h} = \\hat{H}_\\mathrm{tot}^\\ast$という関係で結ばれます. また, [2018年CODATA推奨値](https://physics.nist.gov/cgi-bin/cuu/Value?malsme)より$m_\\mathrm{\\alpha}/m_\\mathrm{e}=7294.299~541~42(24)$が利用できますが, 今回は核の量子効果を無視（$m_\\mathrm{\\alpha}\\rightarrow\\infty$）します. さらに, 以降はこれらのアスタリスクを省略するものとすると, 次のハミルトニアンが得られます.\n", | |
| "\n", | |
| "$$\n", | |
| "\\hat{H} =\n", | |
| "\\hat{H}_0 + \\hat{H}' =\n", | |
| "- \\frac{1}{2} \\nabla_{1}^{2}\n", | |
| "- \\frac{1}{2} \\nabla_{2}^{2}\n", | |
| "- \\frac{2}{r_1}\n", | |
| "- \\frac{2}{r_2}\n", | |
| "+ \\frac{1}{r_{12}}\n", | |
| "$$\n", | |
| "\n", | |
| "このノートでは, これを非摂動ハミルトニアン\n", | |
| "\n", | |
| "$$\n", | |
| "\\hat{H}_0 =\n", | |
| "- \\frac{1}{2} \\nabla_{1}^{2}\n", | |
| "- \\frac{1}{2} \\nabla_{2}^{2}\n", | |
| "- \\frac{2}{r_1}\n", | |
| "- \\frac{2}{r_2}\n", | |
| "$$\n", | |
| "\n", | |
| "と摂動項\n", | |
| "\n", | |
| "$$\n", | |
| "\\hat{H}' =\n", | |
| "\\frac{1}{r_{12}}\n", | |
| "$$\n", | |
| "\n", | |
| "に分けて考えていきます." | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "id": "5574e9dc", | |
| "metadata": {}, | |
| "source": [ | |
| "## 非摂動のエネルギーと波動関数\n", | |
| "\n", | |
| "摂動論では非摂動の問題は既に解けているという前提から始まりますので, まずは非摂動（電子間反発を無視した）ハミルトニアン$\\hat{H}_0$の固有値と固有関数を用意します. よく見ると非摂動ハミルトニアン$\\hat{H}_0$は2つのヘリウムイオンHe⁺のハミルトニアン\n", | |
| "\n", | |
| "$$\n", | |
| "\\hat{h}_1 =\n", | |
| "- \\frac{1}{2} \\nabla_{1}^2\n", | |
| "- \\frac{2}{r_1}\n", | |
| "$$\n", | |
| "\n", | |
| "$$\n", | |
| "\\hat{h}_2 =\n", | |
| "- \\frac{1}{2} \\nabla_{2}^2\n", | |
| "- \\frac{2}{r_2}\n", | |
| "$$\n", | |
| "\n", | |
| "の和になっています. これらの固有値をそれぞれ$E_1$と$E_2$, 固有関数をそれぞれ$\\psi_1(\\pmb{r}_1)$と$\\psi_2(\\pmb{r}_2)$とすると, [2次元の井戸型ポテンシャルや2次元の調和振動子の例](https://www.mns.kyutech.ac.jp/~okamoto/education/quantum/quantum-2dim090603.pdf)などの類推から直ちに非摂動ハミルトニアン$\\hat{H}_0$の固有値は$E_1 + E_2$, 固有関数は$\\psi_1(\\pmb{r}_1) \\psi_2(\\pmb{r}_2)$であると予想されます. 実際, \n", | |
| "\n", | |
| "$$\n", | |
| "\\hat{h}_1 \\psi_1 = E_1 \\psi_1\n", | |
| "$$\n", | |
| "\n", | |
| "$$\n", | |
| "\\hat{h}_2 \\psi_2 = E_2 \\psi_2\n", | |
| "$$\n", | |
| "\n", | |
| "の2式を用いれば\n", | |
| "\n", | |
| "$$\n", | |
| "\\begin{aligned}\n", | |
| "\\hat{H}_0 \\psi_1 \\psi_2\n", | |
| "&=\n", | |
| "\\left( \\hat{h}_1 + \\hat{h}_2 \\right) \\psi_1 \\psi_2 \\\\\n", | |
| "&=\n", | |
| "\\hat{h}_1 \\psi_1 \\psi_2 + \n", | |
| "\\hat{h}_2 \\psi_1 \\psi_2 \\\\\n", | |
| "&=\n", | |
| "\\psi_2 \\hat{h}_1 \\psi_1 + \n", | |
| "\\psi_1 \\hat{h}_2 \\psi_2 \\\\\n", | |
| "&=\n", | |
| "\\psi_2 E_1 \\psi_1 + \n", | |
| "\\psi_1 E_2 \\psi_2 \\\\\n", | |
| "&=\n", | |
| "E_1 \\psi_1 \\psi_2 + \n", | |
| "E_2 \\psi_1 \\psi_2 \\\\\n", | |
| "&=\n", | |
| "\\left( E_1 + E_2 \\right) \\psi_1 \\psi_2\n", | |
| "\\end{aligned}\n", | |
| "$$\n", | |
| "\n", | |
| "となり, 非摂動ハミルトニアン$\\hat{H}_0$について$\\psi_1 \\psi_2$が固有関数となり, $E_1 + E_2$が固有値となることが示せます.\n", | |
| "\n", | |
| "しかし, Schrödinger方程式を満足するだけではダメです. 電子（フェルミ粒子）が2つ以上含まれる系では, それらの空間座標およびスピン座標の置換に対して波動関数は反対称にならなければいけません. これはSchrödinger方程式とは独立した公理であり, Schrödinger方程式を満足するだけは自動的に満足されないため, あらかじめ対称性を満足するような基底関数のみを用いるなどして計算を行わなければいけません. この辺りの話は[A.ザボ, N.S.オストランド著, 大野公男, 阪井健男, 望月祐志訳『新しい量子化学 上 電子構造の理論入門』(東京大学出版会, 1987)](http://www.utp.or.jp/book/b302128.html)の2.1.3に書いてあります. とはいえ, 今回は$n_1=n_2=1$に限定しており, 空間部分が対称でスピン部分が反対称（パラヘリウム）の場合を考えているので$\\psi_1(\\pmb{r}_1) \\psi_2(\\pmb{r}_2)$で問題ありません. 規格化因子が間違っている気がしますが, [Wikipedia](https://ja.wikipedia.org/wiki/%E3%83%98%E3%83%AA%E3%82%A6%E3%83%A0%E5%8E%9F%E5%AD%90#%E9%9B%BB%E5%AD%90%E9%96%93%E5%8F%8D%E7%99%BA%E3%82%92%E7%84%A1%E8%A6%96%E3%81%99%E3%82%8B)に同様の説明が書いてあります.\n", | |
| "\n", | |
| "さて, 具体的なエネルギーを求めましょう. 単純な和ですので, 下記のように計算できます." | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 7, | |
| "id": "59729219", | |
| "metadata": {}, | |
| "outputs": [ | |
| { | |
| "name": "stdout", | |
| "output_type": "stream", | |
| "text": [ | |
| "E = -4.0 Eₕ\n", | |
| " = -108.845544983952 eV\n" | |
| ] | |
| } | |
| ], | |
| "source": [ | |
| "Eₕ2eV = 27.211386245988 # eV # https://physics.nist.gov/cgi-bin/cuu/Value?hrev\n", | |
| "println(\"E = \", (He⁺.E(n=1) + He⁺.E(n=1)), \" Eₕ\")\n", | |
| "println(\" = \", (He⁺.E(n=1) + He⁺.E(n=1)) * Eₕ2eV, \" eV\")" | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "id": "ad558108", | |
| "metadata": {}, | |
| "source": [ | |
| "次に波動関数の規格化条件を確認しましょう. 波動関数は6変数関数として次のように定義します." | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 8, | |
| "id": "5b5e3138", | |
| "metadata": {}, | |
| "outputs": [ | |
| { | |
| "data": { | |
| "text/plain": [ | |
| "Ψ (generic function with 1 method)" | |
| ] | |
| }, | |
| "execution_count": 8, | |
| "metadata": {}, | |
| "output_type": "execute_result" | |
| } | |
| ], | |
| "source": [ | |
| "Ψ(r₁,θ₁,φ₁,r₂,θ₂,φ₂) = He⁺.ψ(r₁,θ₁,φ₁) * He⁺.ψ(r₂,θ₂,φ₂)" | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "id": "60c1e759", | |
| "metadata": {}, | |
| "source": [ | |
| "先ほどと同様にモンテカルロ積分を用いて規格化条件を確認しましょう." | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 9, | |
| "id": "303d4c08", | |
| "metadata": {}, | |
| "outputs": [ | |
| { | |
| "name": "stdout", | |
| "output_type": "stream", | |
| "text": [ | |
| "<ψ|ψ> = 1.0004144350516748 ± 0.0017535335574600287\n" | |
| ] | |
| } | |
| ], | |
| "source": [ | |
| "res = vegas(x ->\n", | |
| " x[1]^2 * sin(x[2]) *\n", | |
| " x[4]^2 * sin(x[5]) *\n", | |
| " abs(Ψ(x...))^2 ,\n", | |
| " [0,0,0,0,0,0],\n", | |
| " [50,π,2π,50,π,2π],\n", | |
| " atol=1e-5\n", | |
| " )\n", | |
| "\n", | |
| "println(\"<ψ|ψ> = \", res.integral_estimate, \" ± \", res.standard_deviation)" | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "id": "0e166092", | |
| "metadata": {}, | |
| "source": [ | |
| "おおよそ1となることが確かめられました." | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "id": "f31afa40", | |
| "metadata": {}, | |
| "source": [ | |
| "## 第一近似のエネルギー\n", | |
| "\n", | |
| "先ほど用意した非摂動波動関数$\\psi_1 \\psi_2$は非摂動ハミルトニアンの厳密な固有関数になっていますが, 摂動項を加えたハミルトニアンの固有関数にはなっていません. しかし, まあ近い形になっているだろうということで$\\psi_1 \\psi_2$を使って摂動項（電子間反発）の期待値を計算します. [こちらの講義ノートの14.3.2](http://www.th.phys.titech.ac.jp/~muto/lectures/QMII09/QMII09_chap14.pdf)に書いてある第一近似のエネルギーの計算に相当します.\n", | |
| "\n", | |
| "波動関数が球面座標で定義されているのでデカルト座標に変換してから電子間距離を計算します.（球面座標から直接計算する方法もあるかもしれません.）" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 10, | |
| "id": "82b91481", | |
| "metadata": {}, | |
| "outputs": [ | |
| { | |
| "data": { | |
| "text/plain": [ | |
| "r₁₂ (generic function with 1 method)" | |
| ] | |
| }, | |
| "execution_count": 10, | |
| "metadata": {}, | |
| "output_type": "execute_result" | |
| } | |
| ], | |
| "source": [ | |
| "x(r,θ,φ) = r*sin(θ)*cos(φ)\n", | |
| "y(r,θ,φ) = r*sin(θ)*sin(φ)\n", | |
| "z(r,θ,φ) = r*cos(θ)\n", | |
| "\n", | |
| "function r₁₂(r₁,θ₁,φ₁,r₂,θ₂,φ₂)\n", | |
| " return sqrt(\n", | |
| " (x(r₁,θ₁,φ₁) - x(r₂,θ₂,φ₂)) ^ 2 + \n", | |
| " (y(r₁,θ₁,φ₁) - y(r₂,θ₂,φ₂)) ^ 2 + \n", | |
| " (z(r₁,θ₁,φ₁) - z(r₂,θ₂,φ₂)) ^ 2\n", | |
| " )\n", | |
| "end" | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "id": "d2a02f42", | |
| "metadata": {}, | |
| "source": [ | |
| "これを先ほどの規格化条件の確認をした積分のコードに挟み込むと, 摂動項（電子間反発）の期待値が計算できます." | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 11, | |
| "id": "0f5a0cde", | |
| "metadata": {}, | |
| "outputs": [ | |
| { | |
| "name": "stdout", | |
| "output_type": "stream", | |
| "text": [ | |
| "<V₁₂> = 1.2495678829756485 Eₕ\n", | |
| " = 34.0024743042319 eV\n" | |
| ] | |
| } | |
| ], | |
| "source": [ | |
| "res = vegas(x ->\n", | |
| " x[1]^2 * sin(x[2]) *\n", | |
| " x[4]^2 * sin(x[5]) *\n", | |
| " conj(Ψ(x...)) *\n", | |
| " 1 / r₁₂(x...) *\n", | |
| " Ψ(x...) ,\n", | |
| " [0,0,0,0,0,0],\n", | |
| " [50,π,2π,50,π,2π],\n", | |
| " rtol = 1e-6,\n", | |
| " atol = 1e-6\n", | |
| " )\n", | |
| "\n", | |
| "∫ψψV₁₂ψψdr₁dr₂ = res.integral_estimate\n", | |
| "\n", | |
| "println(\"<V₁₂> = \", ∫ψψV₁₂ψψdr₁dr₂, \" Eₕ\")\n", | |
| "println(\" = \", ∫ψψV₁₂ψψdr₁dr₂ * Eₕ2eV, \" eV\")" | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "id": "b54df961", | |
| "metadata": {}, | |
| "source": [ | |
| "これを先ほどの$-4 E_\\mathrm{h}$に足すと, 摂動項も含めた全エネルギーが次のように計算できます." | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 12, | |
| "id": "3775a8fd", | |
| "metadata": {}, | |
| "outputs": [ | |
| { | |
| "name": "stdout", | |
| "output_type": "stream", | |
| "text": [ | |
| "<H> = -2.7504321170243515 Eₕ\n", | |
| " = -74.8430706797201 eV\n" | |
| ] | |
| } | |
| ], | |
| "source": [ | |
| "E = He⁺.E(n=1) + He⁺.E(n=1) + ∫ψψV₁₂ψψdr₁dr₂\n", | |
| "\n", | |
| "println(\"<H> = \", E, \" Eₕ\")\n", | |
| "println(\" = \", E * Eₕ2eV, \" eV\")" | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "id": "da4877a5", | |
| "metadata": {}, | |
| "source": [ | |
| "これは現在知られている最も正確なエネルギー [$-2.903~724~377 034~119~598~311~159~245~194~4~E_\\mathrm{h}$](https://doi.org/10.1103/PhysRevLett.101.240406)と比べると高いですが, 全くおかしな値ではないことが確認できます.\n", | |
| "\n", | |
| "変分原理の観点から言えば, 摂動項も含めたハミルトニアンにとって非摂動波動関数は的外れではないものの, あまり正確な波動関数にはなっていないであろうということです. なお, 変分原理の意味するところは「エネルギーが低くければ低くいほど**エネルギーが正確**である」という意味です. 他の期待値も正確になるとは限りませんので, 波動関数の良し悪しを判断する絶対的な基準ではないということに注意してください. 詳しくは[こちらの解説](https://doi.org/10.1021/acs.jchemed.8b00959)をご覧ください." | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "id": "747049c3", | |
| "metadata": {}, | |
| "source": [ | |
| "## おわりに\n", | |
| "\n", | |
| "摂動論は第一近似だけ考えて終わりではなく, より高次の補正を加えていくアプローチです. このノートの内容だけだとHartree積を試行波動関数に用いて変分計算を行ったこととの区別がつかないため, どちらを強調するか迷いました. ここでは非摂動波動関数を用いて摂動項の期待値を計算することの正当性が確認できたということで満足するものとします." | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "id": "55868543", | |
| "metadata": {}, | |
| "source": [ | |
| "## 参考文献\n", | |
| "\n", | |
| "- [東京工業大学の武藤先生の講義ノート](http://www.th.phys.titech.ac.jp/~muto/lectures/QMII09/QMII09_chap14.pdf)\n", | |
| "- [Antique.jl](https://ohno.github.io/Antique.jl/dev/)\n", | |
| "- [MonteCarloIntegration.jl](https://github.com/ranjanan/MonteCarloIntegration.jl)\n", | |
| "- [IUPAC Greenbook 3.2.9](https://unit.aist.go.jp/nmij/public/report/others/pdf/iupac_green_book_jp.pdf)\n", | |
| "- [2018年CODATA推奨値](https://physics.nist.gov/cgi-bin/cuu/Value?hrev)\n", | |
| "- [A.ザボ, N.S.オストランド著, 大野公男, 阪井健男, 望月祐志訳『新しい量子化学 上 電子構造の理論入門』(東京大学出版会, 1987)](http://www.utp.or.jp/book/b302128.html)\n", | |
| "- [H. Nakashima, H. Nakatsuji, _Phys. Rev. Lett._ **101**, 240406 (2008)](https://doi.org/10.1103/PhysRevLett.101.240406)\n", | |
| "- [J. F. Pérez-Torres, _J. Chem. Educ._ **96**, 4, 704–707 (2019)](https://doi.org/10.1021/acs.jchemed.8b00959)" | |
| ] | |
| } | |
| ], | |
| "metadata": { | |
| "kernelspec": { | |
| "display_name": "Julia 1.9.2", | |
| "language": "julia", | |
| "name": "julia-1.9" | |
| }, | |
| "language_info": { | |
| "file_extension": ".jl", | |
| "mimetype": "application/julia", | |
| "name": "julia", | |
| "version": "1.9.2" | |
| } | |
| }, | |
| "nbformat": 4, | |
| "nbformat_minor": 5 | |
| } |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment