李媛媛,張慶偉
(1.北京信息科技大學 理學院,北京 100192;2.北京奇虎360科技有限公司,北京 100015)
狀態(tài)方程[1]是描述物質系統(tǒng)中各狀態(tài)變量之間關系的函數(shù)表達式,用來表達在一定熱力學條件下物質的性質。通常人們研究最多的是力學響應特性的狀態(tài)方程,即所謂的狹義狀態(tài)方程。溫稠密物質[2-3]是一種不常見物質,廣泛存在于地幔內部和宇宙星體內,具有高的能量密度特性。溫稠密物質中存在復雜的電離效應,精確了解不同粒子的電離程度,可以更好地了解強耦合下溫稠密物質內各種粒子和束團的狀態(tài)與成分。
在科學研究中,想要獲得溫稠密物質的狀態(tài)方程十分困難。只有在超高壓超高溫的情況下,才會出現(xiàn)溫稠密物質,一般需要通過核聚變實驗進行相關研究,但是花費巨大,所以目前許多研究學者已經使用遺傳算法來獲得狀態(tài)方程。王成[4]提出基于基因遺傳算法和γ律狀態(tài)方程的JWL狀態(tài)方程參數(shù)的確定;溫麗晶[5]采用遺傳算法確定PBXC03等炸藥轟爆產物JWL狀態(tài)方程的參數(shù),應用確定的狀態(tài)方程參數(shù)模擬爆轟驅動遠觀過程,計算結果與實驗結果符合。這2篇文獻的算法都沒有采用并行,效率較低,精確度較差。張慶偉[6]提出的基于Spark的并行遺傳算法得到了較為理想的效果。本文使用基于Spark并行遺傳算法來確定溫稠密物質狀態(tài)方程的參數(shù)。
本次實驗的數(shù)據(jù)選自王聰?shù)冉⒌臍?氘、氚)的密度泛函第一性原理寬區(qū)狀態(tài)方程數(shù)據(jù)庫[7]。數(shù)據(jù)主要為氫、氘、氚在極端條件下的壓強、密度、溫度和內能。其中的密度范圍為10-3~103g/cm3,溫度范圍為103~107K。按照密度劃分為31組,每組有13~21條數(shù)據(jù)。本文的研究都是基于表1的數(shù)據(jù)。
表1 氫在極端條件下的數(shù)據(jù)
在自然界中很多物質的狀態(tài)方程是未知的,所以本文需要先根據(jù)實驗數(shù)據(jù)確定狀態(tài)方程的形式。在密度固定時,壓力和溫度之間存在著一定的函數(shù)關系,假設兩者之間的函數(shù)關系為
f(T,P)=0
(1)
表1中每組數(shù)據(jù)為1321條不等,而泰勒的三階展開需要10個參數(shù),四階展開需要15個參數(shù),為了可以求解到一組唯一解,選擇三階展開,式(1)在(0,0)處的三階泰勒展開式為
f(T,P)=A+BT+CP+DT2+ETP+FP2+
GT3+HT2P+ITP2+JP3+ο3
(2)
式中:A,B,C,D,E,F,G,H,I,J為狀態(tài)方程對應的參數(shù);T、P分別為溫度、壓強。
考慮到隨著密度的增大,壓力增大的速度極快,可以達到107,為了方便計算,對T和P分別取對數(shù),令t=log10T,p=log10P,代入式(2),得
f(t,p)=A+Bt+Cp+Dt2+Etp+Fp2+
Gt3+Ht2p+Itp2+Jp3+ο3
(3)
從31組不同密度的數(shù)據(jù)中,每組選10個數(shù)據(jù)代入式(3),每10組數(shù)據(jù)都會得出一個式(4),再對式(4)求解非線性方程組即可以求出一組唯一的參數(shù)(A,B,C,D,E,F,G,H,I,J)。
f(t,p)=
(4)
采用矩陣法可以求出對應31組密度的溫稠密物質狀態(tài)方程參數(shù)。表2為部分結果。
表2顯示不同的密度對應的溫稠密物質狀態(tài)方程的參數(shù),考慮是否可以尋求到ρ和每個參數(shù)之間的函數(shù)關系。令:
(5)
式中(ma0,ma1,ma2,ma3,ma4)為參數(shù)A在不同密度下所對應的的參數(shù)。同理,對于所有的參數(shù),都有一組(m0,m1,m2,m3,m4) 與之對應。其矩陣式為
(6)
式中[An,Bn,Cn,Dn,En,Fn,Gn,Hn,In,Jn] (n∈[1,31])對應著表2中節(jié)選的部分數(shù)據(jù)。可求出每個參數(shù)對應的一組m值,如表3所示。表3為使用ρ的對數(shù)進行計算得出溫稠密物質狀態(tài)方程參數(shù)在不同密度下所對應的參數(shù)。把表2中得到的狀態(tài)方程的參數(shù)代入式(7)中即可求出溫稠密物質狀態(tài)方程:
A+Bt+Cp+Dt2+Etp+Fp2+Gt3+
Ht2p+Itp2+Jp3+ο3=0
(7)
表2 矩陣法求得狀態(tài)方程參數(shù)的部分結果
表3 參數(shù)結果矩陣法求解的狀態(tài)方程參數(shù)在不同密度下對應的參數(shù)
矩陣法可以求出相應的狀態(tài)方程,但是存在一定的缺陷。因為不同密度對應的數(shù)據(jù)量不同,上面的算法只選用了其中的10組數(shù)據(jù),雖然可以求出一組確定解,但是該解對未參與計算的數(shù)據(jù)的適用性有一定的不確定性。雖然可以通過后期驗證其適用性,但是解一旦確定后就無法修改。如果求得的解不適應未參與計算的數(shù)據(jù),則會導致該解無效。所以針對上述方法的缺陷,本文使用基于Spark的并行遺傳算法來確定狀態(tài)方程的參數(shù)。
張慶偉提出基于Spark的并行遺傳算法可以解決理想氣體的狀態(tài)方程的參數(shù),并且計算結果良好,在計算時間和結果的精度上都取得了較好的效果。本文將該方法用于確定溫稠密物質狀態(tài)方程參數(shù),使用Scala[8-9]語言編寫基于Spark的并行遺傳算法程序,將其在一臺服務器的虛擬機上運行。Leque[10]介紹了3種遺傳算法的并行方式:分布式模型、蜂窩煤模型和混合模型。Spark是一種主從式的分布式計算框架,適用于全局分布式的并行遺傳算法。但傳統(tǒng)的遺傳算法不能直接結合spark來確定狀態(tài)方程的參數(shù),需要對傳統(tǒng)遺傳算法進行改進。主要針對以下幾個方面。
1)傳統(tǒng)遺傳算法始終在一個固定的空間內進行迭代,算法效率低下,基于spark并行遺傳算法將固定的空間設計為變化的搜索空間。具體方法為:將最初的搜索空間劃分為若干個子區(qū)間,且每個子區(qū)間對應不同的子種群,將子種群及其對應子空間放到不同的節(jié)點上進行運算,選擇得到最小適應度個體所在的子空間,再將其繼續(xù)劃分,達到不斷變化搜索空間的效果。
2)因為基于spark并行遺傳算法需要不斷地變化搜索空間,所以變化的子空間的子種群也需要同步更新,然后再對狀態(tài)方程未知參數(shù)的取值范圍進行劃分,在每個子空間上隨機生成子種群,再將生成好的種群寫入文件中,每行代表一個區(qū)間的所有個體及其所在的取值空間,以便后續(xù)使用,此時子種群已生成。
3)在產生子種群后生成新種群,然后對新種群個體進行適應度評價,通過適應度函數(shù)方法得到新種群的最小適應度、對應區(qū)間、最小適應度對應的個體,再把所有的數(shù)據(jù)回收到主節(jié)點,在主節(jié)點上選擇最小適應度值所在的取值空間,將空間重新劃分,產生子種群,將產生的子種群繼續(xù)劃分若干子種群,直至滿足截止條件,將適應度最小的個體保存下來,并將其加入新的子種群中,這樣就重新選定了搜索空間。
4)傳統(tǒng)遺傳算法有的只并行適應度函數(shù)的計算,有的只并行遺傳操作。為了將遺傳算法和spark相結合來確定狀態(tài)方程中的參數(shù),將傳統(tǒng)的串行遺傳算法[11-12]進行了改進:①算法中包含2個函數(shù)并行運行,一個是遺傳操作的并行,另一個是獲得個體適應度函數(shù)的并行運算。因為產生新種群替代舊種群需要的循環(huán)次數(shù)多,且選擇、交叉、變異等操作也比較費時,所以遺傳比獲得個體適應度需要更多的時間成本;②針對遺傳算法中的迭代操作。首先初始化迭代次數(shù),判斷是否滿足截止條件,滿足返回新種群,否則進入迭代,從Pstring[13-14]中獲得狀態(tài)方程的參數(shù)區(qū)間,獲得種群,判斷新種群個體數(shù)是否等于舊種群個體數(shù),如果達到就將舊種群中適應度最小的個體替代新種群中適應度最大的個體,進行下一次迭代;否則就選出2個個體,進行遺傳操作,將新產生的個體加入到新種群中,直到新種群中的個體數(shù)滿足條件。在基于Spark的并行遺傳算法中,因為狀態(tài)方程形式未知,所以遺傳算法中適應度函數(shù)的確定難度變大,根據(jù)表1中數(shù)據(jù),由于溫度和壓強存在一定的函數(shù)關系,令
(8)
則適應度函數(shù)可以寫為
(9)
式中n為每組數(shù)據(jù)的數(shù)量。
表1中的t、p值是通過物理實驗測得,將t值和用遺傳算法產生的p值代入式(9),根據(jù)適應度函數(shù)的變化情況確定狀態(tài)方程的參數(shù)。表5為所求狀態(tài)方程參數(shù)的部分結果。
表5 遺傳算法所得狀態(tài)方程參數(shù)的部分結果
與矩陣法相同,遺傳算法對于所求得的參數(shù)和密度之間也存在一定的函數(shù)關系,令
(10)
式中γ1、γ2、γ3、γ4、γ5分別為參數(shù)A、B、C、D、E在不同密度下所對應的參數(shù)。
適應度函數(shù)為
(11)
式中z為數(shù)據(jù)分組數(shù)。使用遺傳算法求出參數(shù)在不同密度下所對應的參數(shù),如表6所示。
表6 遺傳算法求解的狀態(tài)方程參數(shù)在不同密度下對應的參數(shù)
本文中的矩陣法和基于Spark的并行遺傳算法使用的都是表1中的數(shù)據(jù)。2種方法在不同密度下誤差的變化趨勢如圖1所示。
圖1 不同分組的實驗結果誤差
從圖1可以看出,使用矩陣法所得的參數(shù)較為精確,使用遺傳算法所得結果的精度較矩陣法有所下降,且隨著密度的增大,遺傳算法求出的數(shù)據(jù)與物理實驗測得數(shù)據(jù)的誤差也越來越大。主要原因是隨著密度的增加,壓強變化范圍也逐漸增加,且增長速度很快,過大的范圍容易產生狀態(tài)方程的分區(qū)域問題,即不同的溫度和壓強區(qū)域狀態(tài)方程的參數(shù)和結構都會發(fā)生變化,導致誤差逐漸增大。
圖2是2種方法在相同密度下所求參數(shù)的平均誤差。從圖2可以看出,矩陣法對應的10個參數(shù)的平均誤差不夠穩(wěn)定,尤其是A、B、C誤差過大,根本無法滿足要求。而遺傳算法對應參數(shù)的平均誤差極小且十分穩(wěn)定。
圖2 不同參數(shù)對應的誤差
結合圖1和圖2可以看出,對于不同密度分組下的參數(shù)求解,矩陣法求得溫稠密物質狀態(tài)方程參數(shù)的誤差在密度偏大的時候要優(yōu)于遺傳算法的結果;但是,對于在相同密度時,求溫稠密物質狀態(tài)方程參數(shù)平均誤差時,遺傳算法在精度和穩(wěn)定性上都要優(yōu)于矩陣法。此外,矩陣法所求的參數(shù)對數(shù)據(jù)的依賴性更大,只適用某些數(shù)據(jù),該解對未參與計算的數(shù)據(jù)適用性有一定的不確定性;而遺傳算法的解基本可以忽略與物理實驗的數(shù)據(jù)誤差,更加精確。
為了盡可能減少實驗次數(shù),利用更少的數(shù)據(jù)獲得更加準確的狀態(tài)方程,本文采用基于Spark的并行遺傳算法構建溫稠密物質狀態(tài)方程。實驗結果表明,本文方法所求溫稠密物質狀態(tài)方程參數(shù)的穩(wěn)定性優(yōu)于矩陣法。對于密度逐漸增大時,遺傳算法會有較大誤差的問題,可以考慮把不同溫度、密度和壓強進行分區(qū)域研究,以便得到更精確的結果。