深度學習與量子蒙特卡洛 (3) 最優路徑與神經網路

05-18

來自專欄量子蒙特卡洛與機器學習

0) 讓我們對之前文章的討論做一個小小的總結。記住我們的目標是計算配分函數

$Z = int_{-infty}^{infty} dx~ ext{e}^{-S}.~~~~~~~~~~~~~~~~~( ext{Eq.}~1)$

當作用量 $S$ 為實數時，該積分可以通過Monte-Carlo方法計算出來；當 $S$ 的虛部不為零時Monte-Carlo方法會失效。為了度量失效的程度，我們引入積分：

$Z_ ext{ph} =int_{-infty}^{infty} dx~ ext{e}^{-S_{Re}}= int_{-infty}^{infty} dx~ ext{e}^{-Re(S)}~~~~~~~~( ext{Eq}. 2)$

記 $Z_ ext{ph}$ 為受限配分函數。定義受限因子

$heta = Z/Z_ ext{ph}~~~~~~~~( ext{Eq}. 3)$

$heta$ 度量了作用量虛部的效應。當 $heta<10^{-2}$ 時，我們稱積分(Eq.1)存在符號問題；當 $heta<10^{-4}$ 時，我們稱積分(Eq.1)存在嚴重的符號問題；當 $heta<10^{-8}$ 時，我們稱積分(Eq.1)存在極其嚴重的符號問題；

1) 當作用量 $S$ 的虛部不為零時( $Im(S) eq 0$ ), 積分(Eq.1)無法用Monte-Carlo方法計算，但積分(Eq.2)可以用Monte-Carlo方法計算。

$heta = frac{Z}{Z_ ext{ph}} = frac{int dx~e^{-i~Im(S)}e^{-Re(S)}}{int dx~e^{-Re(S)} } = <e^{-i~Im(S)}>_ ext{ph}~~~~~~~~( ext{Eq}. 4)$

任意物理量 $mathcal{O}$ 的期望值可以用如下公式計算：

$<mathcal{O}> = frac{int dx~mathcal{O}e^{-i~Im(S)}e^{-Re(S)}}{int dx~e^{-S} } = frac{<mathcal{O}e^{-i~Im(S)}>_ ext{ph}}{ heta} ~~~~~~~~( ext{Eq}. 5)$

2) 考慮如下作用量

$S=x^2-p~ ext{log}(x/a+i)~~~~~~~~~~~~~~~~~~~~~( ext{Eq}. 6)$

將積分路徑延拓到複平面，這裡我們只做一個簡單的平移 $x ightarrow z = x + ieta ~(x,etain R)$

利用公式（Eq.4）計算受限因子 $heta$ 隨著 $eta$ 的變化，得到下圖：

可以看出當 $heta<10^{-2}$ 時蒙特卡洛方法給出的值開始和真實值發生偏離。計算中取了300000個點，如果增加投點數目，精度是可以增加的。注意在上圖中我們計算的是 $heta$ 的模。

如何使用神經網路尋找最優參數？待續（18/01/09 11：45）

附錄

python code:

# written by K. J. Sun # copy right reservedfrom __future__ import print_functionimport numpy as npimport matplotlib.pyplot as pltimport matplotlib.image as mpingfrom pyhmc import hmcdef logprob(t, beta): alpha = 10 p = 50 x = t+ 1j*beta S = (x)**2/2.-p*np.log(x+1j*alpha)+p*np.log(alpha) logp=-np.real(S) grad = -np.real(2*x-p/(x+1j*alpha)) return logp, grad def toy(O,t, beta): alpha = 10 p = 50 x = t+ 1j*beta S = (x)**2/2.-p*np.log(x+1j*alpha)+p*np.log(alpha) f = (x+1j*alpha)**p*np.exp(-x**2/2)/alpha**p phase = f/abs(f) Ob = phase*O #fabs = abs(f) return phase,Ob betav = np.linspace(-10.0,0,50)k = 0phasev = np.zeros(50)for beta in betav:#beta = -5. print("step= "+str(k+1)) samples1 = hmc(logprob, x0=np.array([5.0]), args=(beta,), n_samples=300000) samples1 = samples1[10000:-1] samples2 = -samples1 #hmc(logprob, x0=np.array([-5.0]), args=(beta,), n_samples=100000) samples = np.row_stack((samples1,samples2)) phase,Ob = toy((samples+1j*beta)**2,samples,beta)#Zabs = np.mean(fabs) avphase = np.mean(phase) O = np.mean(Ob)/avphase#Z = Zabs*avphase#print(Zabs) print(abs(avphase))#print(np.real(Z)) print(np.real(O)) phasev[k] = abs(avphase) k=k+1plt.figure()plt.plot(betav,phasev)plt.show()

save it as hmc_toy.py and run it python3 hmc_toy.py. It takes about 30*60 seconds.