李國(guó)慶, 孫 威, 王 穎, 劉福窯, 鄭晶晶

(上海工程技術(shù)大學(xué) 數(shù)理與統(tǒng)計(jì)學(xué)院, 上海 201620)

尋求引起全球氣候及地表環(huán)境周期性變化成因的研究一直是熱點(diǎn).行星系運(yùn)動(dòng)作為重要的天文因素,對(duì)于研究全球氣候及地表環(huán)境周期性變化具有重要意義.目前學(xué)術(shù)界普遍接受的觀點(diǎn)認(rèn)為,太陽(yáng)活動(dòng)具有多尺度的周期特征,除普遍發(fā)現(xiàn)的約11 a最重要的Schwabe周期波動(dòng)成分外,還有約22 a的Hale周期、約210 a的Suess周期等不同尺度的周期結(jié)果[1-3],且研究逐步聚焦在太陽(yáng)活動(dòng)周期的觀察、檢測(cè)和預(yù)報(bào).半個(gè)多世紀(jì)以來(lái),關(guān)于太陽(yáng)活動(dòng)及行星運(yùn)動(dòng)的研究取得了頗為豐富的成果.相關(guān)研究成果表明,太陽(yáng)活動(dòng)的諸多表征指標(biāo)與太陽(yáng)黑子具有良好的相關(guān)性,太陽(yáng)黑子活動(dòng)具有準(zhǔn)11 a周期性規(guī)律并且與全球氣候變化存在顯著響應(yīng)[4-7]:一方面,利用長(zhǎng)時(shí)空尺度的太陽(yáng)黑子相對(duì)數(shù)及太陽(yáng)黑子群磁極性更迭研究太陽(yáng)活動(dòng)及太陽(yáng)運(yùn)動(dòng)的周期規(guī)律,不失為行之有效的得力工具[8];另一方面,借助太陽(yáng)磁場(chǎng)22 a的磁性周期,可研究探索太陽(yáng)運(yùn)動(dòng)軌道的周期性規(guī)律[9].劉復(fù)剛等[10]通過(guò)質(zhì)心定義,充分考慮各行星質(zhì)量差異對(duì)行星系統(tǒng)運(yùn)動(dòng)軌跡的影響,確立行星會(huì)合指數(shù)運(yùn)動(dòng)指標(biāo)并給出運(yùn)動(dòng)學(xué)方程,將較為繁雜的天體間繞轉(zhuǎn)模型進(jìn)行簡(jiǎn)化,在一定程度上彌補(bǔ)了僅能夠研究極少數(shù)行星會(huì)合或離散的現(xiàn)狀.孫威等[11-12]考慮到僅依賴行星會(huì)合指數(shù)無(wú)法對(duì)行星系質(zhì)心及太陽(yáng)系質(zhì)心準(zhǔn)確定位,為解決這一問(wèn)題,創(chuàng)建行星系日心經(jīng)度,數(shù)值模擬結(jié)果顯示,天體系統(tǒng)內(nèi)互相關(guān)聯(lián)的繞轉(zhuǎn)系質(zhì)心的運(yùn)轉(zhuǎn)軌跡具有11.5 a的周期性規(guī)律,太陽(yáng)系質(zhì)心相對(duì)于太陽(yáng)本體范圍內(nèi)出現(xiàn)的概率約為34.5%.研究表明,行星系統(tǒng)周期性運(yùn)動(dòng)與太陽(yáng)運(yùn)動(dòng)及太陽(yáng)活動(dòng)的變化規(guī)律存在顯著一致性,太陽(yáng)周期性活動(dòng)勢(shì)必引起全球氣候及地表環(huán)境變化[13].行星系統(tǒng)運(yùn)動(dòng)的周期性規(guī)律對(duì)于表征全球氣候及地表環(huán)境變化的周期性至關(guān)重要.因此,構(gòu)建行星會(huì)合指數(shù),可為研究太陽(yáng)活動(dòng)和太陽(yáng)運(yùn)動(dòng)周期性規(guī)律及其引發(fā)地表環(huán)境及全球氣候變化的周期性規(guī)律提供新的研究思路和途徑.本文借助已構(gòu)建的行星會(huì)合指數(shù)運(yùn)動(dòng)學(xué)方程,利用索周法對(duì)1000—2000AD的基礎(chǔ)數(shù)據(jù)進(jìn)行周期性分析,進(jìn)而研究行星系統(tǒng)運(yùn)動(dòng)特征及周期性變化規(guī)律.

1 行星會(huì)合指數(shù)的建立及數(shù)值模擬

1.1 行星會(huì)合指數(shù)構(gòu)建

根據(jù)靜力矩平衡原理對(duì)質(zhì)心的定義,設(shè)由n個(gè)質(zhì)點(diǎn)組成的質(zhì)量為M的質(zhì)點(diǎn)系,其系統(tǒng)內(nèi)各質(zhì)點(diǎn)的質(zhì)量分別為mp1,mp2,…,mpn.若rp1,rp2,…,rpn分別為質(zhì)點(diǎn)P關(guān)于相對(duì)于固定點(diǎn)(坐標(biāo)原點(diǎn))的矢徑,則質(zhì)點(diǎn)系質(zhì)心C的位矢rC為[14]

(1)

劉復(fù)剛等[10]根據(jù)理論力學(xué)中質(zhì)心方程的定義,類比定義出行星系質(zhì)心為太陽(yáng)系各繞轉(zhuǎn)運(yùn)動(dòng)的行星以各自權(quán)重質(zhì)量為系數(shù),并以其距繞轉(zhuǎn)中心天體的位矢半徑為矢量和.太陽(yáng)系各繞轉(zhuǎn)天體的天文參數(shù)數(shù)據(jù)見表1.表中數(shù)據(jù)來(lái)源于美國(guó)國(guó)家航空航天局(National Aeronautics and Space Administration,NASA).

表1 太陽(yáng)系行星軌道參數(shù)Table 1 Orbit parameters of planet in solar system

利用表1數(shù)據(jù),將平均繞轉(zhuǎn)半徑r轉(zhuǎn)換成以質(zhì)量作權(quán)重系數(shù)的位矢半徑ri,創(chuàng)建行星系質(zhì)心繞日會(huì)合指數(shù)(K)的數(shù)學(xué)表達(dá)式[7]為

K=

(2)

式中:K為行星會(huì)合指數(shù),K≥0;ωi為行星繞轉(zhuǎn)的平均角速度;φi為行星初始日心經(jīng)度;ri為各繞轉(zhuǎn)行星以其質(zhì)量為權(quán)重的位矢半徑;t為時(shí)間.

1.2 行星會(huì)合指數(shù)數(shù)值模擬

通過(guò)NASA平臺(tái)數(shù)據(jù)整理獲得2000年1月1日12時(shí)日心經(jīng)度并作數(shù)據(jù)處理,見表2.結(jié)合表1中各行星軌道天文參數(shù)處理后代入式(2),得到行星會(huì)合指數(shù)變化,如圖1所示.

表2 八大行星的日心經(jīng)度Table 2 Heliocentric longitude of eight planets

圖1 1000—2000AD行星會(huì)合指數(shù)圖Fig.1 Diagram of planet juncture index from 1000 to 2000AD

圖中刻畫出行星系質(zhì)心距太陽(yáng)質(zhì)心的距離大小,顯示的數(shù)值變化曲線可理解為行星系質(zhì)心繞不動(dòng)點(diǎn)的日心軌跡,亦可視為任意“兩心”之間距離的動(dòng)態(tài)變化.經(jīng)理論計(jì)算發(fā)現(xiàn),K的取值為0～7.520 2,但實(shí)際情況中均無(wú)法達(dá)到最大值和最小值.理論上當(dāng)太陽(yáng)系各行星日心經(jīng)度近似相等時(shí),行星會(huì)合指數(shù)達(dá)到最大值;當(dāng)權(quán)重質(zhì)量最大的木星(權(quán)重占比49.151 569%)與其他行星之間相差半圓周日心經(jīng)度時(shí),行星會(huì)合指數(shù)達(dá)到最小值[10].K取最大值和最小值時(shí)分別對(duì)應(yīng)其高度會(huì)合和分離狀態(tài).然而,根據(jù)各行星權(quán)重半徑構(gòu)建的行星會(huì)合指數(shù)運(yùn)動(dòng)學(xué)方程,在某種意義上反映的是地外的木星、土星、天王星、海王星等4顆大質(zhì)量行星會(huì)合程度[10].

2 周期性分析方法及應(yīng)用

2.1 索周法

索周法是通過(guò)數(shù)據(jù)處理、分析其變化特征而獲得相關(guān)數(shù)據(jù)變化規(guī)律的一種周期性分析方法.該方法處理的基本流程為:按照不同類型矩陣排列數(shù)據(jù)并對(duì)矩陣中的各列數(shù)值求和,利用各列數(shù)值之差的絕對(duì)值大小反映數(shù)據(jù)變化的顯隱性.將蘊(yùn)含周期性規(guī)律的數(shù)據(jù)排列在不同行數(shù)或列數(shù)的系列矩陣中,對(duì)各矩陣不同列數(shù)值求和,各求和數(shù)值之差的絕對(duì)值反映在以矩陣列數(shù)代表周期的曲線中,即曲線中的波峰或波谷坐標(biāo)值[15].作為應(yīng)用型周期分析方法,索周法處理連續(xù)時(shí)間序列呈現(xiàn)出普適性、精確性等特點(diǎn),可規(guī)避典型周期處理方式涉及信號(hào)的連續(xù)性、基函數(shù)選擇、分解方式選擇的不確定性,較為全面準(zhǔn)確地分析出不同長(zhǎng)度數(shù)據(jù)所隱含的周期.陳彪等[16]借助1750年后210 a的太陽(yáng)黑子相對(duì)數(shù)月平滑平均值,利用索周法推導(dǎo)求出太陽(yáng)黑子相對(duì)數(shù)曲線的主、副周期[16],根據(jù)圖像顯示的周期特點(diǎn)成功建構(gòu)出相對(duì)數(shù)曲線的近似表達(dá)式.本文借助行星會(huì)合指數(shù)模擬數(shù)據(jù),結(jié)合索周法探求其變化的周期性韻律,考慮索周矩陣行數(shù)和列數(shù)及數(shù)據(jù)排列方式等,為使獲取的周期結(jié)果更為準(zhǔn)確,確定索周法的具體步驟如下.

1) 將行星會(huì)合指數(shù)方程數(shù)值模擬獲得的基礎(chǔ)數(shù)據(jù),從1500年開始按照步長(zhǎng)為1 a取K值.

2) 將會(huì)合指數(shù)數(shù)值按照相同時(shí)間間隔依次排列成列數(shù)為8列,9列,…,50列,行數(shù)為10行,11行,12行的一一對(duì)應(yīng)的矩陣.行星會(huì)合指數(shù)在1500—1600AD,按照10行10列矩陣排列,即T=10,行星會(huì)合指數(shù)排列矩陣見表3.

3) 選取行數(shù)相同列數(shù)不同的矩陣,將不同矩陣的每列數(shù)值相加得最大值和最小值,稱這兩個(gè)最值數(shù)值差的絕對(duì)值為Δ(表3中,Δ=43.310 4-37.764 5=5.545 9).

表3 行星會(huì)合指數(shù)排列矩陣Table 3 Numerical arrangement matrix of planet juncture index

4) 以T為橫坐標(biāo),Δ為縱坐標(biāo)繪制曲線(正排、逆排).利用歸一化數(shù)據(jù)處理方式,將不同列數(shù)矩陣所對(duì)應(yīng)的Δ值反映在以矩陣列數(shù)數(shù)值代表數(shù)據(jù)隱含周期的坐標(biāo)系中,即列數(shù)數(shù)值代表周期數(shù)值,其中曲線峰值(Δ值)代表周期性結(jié)果的強(qiáng)弱.

2.2 行星會(huì)合指數(shù)的周期性分析

太陽(yáng)周期性活動(dòng)在調(diào)制北大西洋濤動(dòng)、響應(yīng)北半球平流層溫度異常和驅(qū)動(dòng)地球氣候等影響氣候變化和地表環(huán)境研究中日益受到關(guān)注[17-19].對(duì)已研究發(fā)現(xiàn)的太陽(yáng)軌道周期性規(guī)律與表征環(huán)境變化現(xiàn)象特征進(jìn)行關(guān)聯(lián)研究,嘗試將行星運(yùn)動(dòng)與驅(qū)動(dòng)全球氣候變化的內(nèi)因組建關(guān)聯(lián).通過(guò)K運(yùn)動(dòng)學(xué)表達(dá)式,能夠巧妙地運(yùn)用太陽(yáng)系各行星運(yùn)動(dòng)的變化規(guī)律,轉(zhuǎn)換為探討系統(tǒng)內(nèi)行星系質(zhì)心繞太陽(yáng)運(yùn)動(dòng)的周期性規(guī)律;同時(shí),行星會(huì)合指數(shù)運(yùn)動(dòng)學(xué)方程的數(shù)值模擬結(jié)果顯示,太陽(yáng)繞太陽(yáng)系質(zhì)心運(yùn)動(dòng)由19.858 5 a的基本周期構(gòu)成[10].

根據(jù)行星會(huì)合指數(shù)數(shù)值模擬時(shí)間序列,合理劃分涵蓋不同峰值、谷值的周期區(qū)段,如圖2所示.截取3段不同區(qū)間長(zhǎng)度的整周期片段.其中,TM6-TM5=179.21 a、TM4-TM3=180.30 a、TM2-TM1=179.35 a即平均周期為(179.21+180.30+179.35)/3=179.62 a.與此同時(shí),該整周期片段包含9個(gè)小周期,每個(gè)小周期長(zhǎng)度為19.957 7 a,該數(shù)值與太陽(yáng)系兩質(zhì)心互相繞轉(zhuǎn)周期19.858 5 a大致吻合[10].

圖2 行星系質(zhì)心運(yùn)動(dòng)特征及其周期Fig.2 Movement characteristics of mass center and its period of planetary system

利用表1和表2數(shù)值進(jìn)行模擬,獲取1000—1500AD和1500—2000AD上K值.依照索周步驟對(duì)基礎(chǔ)數(shù)據(jù)依次處理,考慮矩陣排列方式及尾排數(shù)據(jù)對(duì)結(jié)果精準(zhǔn)度的影響,選取行數(shù)為10行(定值)、列數(shù)不同的矩陣,采取正、逆排列數(shù)據(jù)方式進(jìn)行處理.將行數(shù)均為10行、列數(shù)為8～50列的矩陣正逆排列,對(duì)各矩陣列數(shù)值求和并求最大值與最小值差值的絕對(duì)值Δ,將不同列數(shù)矩陣所得Δ數(shù)據(jù)歸一化處理(Δ/Δ0),反映在以不同列數(shù)數(shù)值代表數(shù)據(jù)隱含周期的坐標(biāo)系中,如圖3所示.圖中,曲線峰值代表周期性規(guī)律的強(qiáng)弱程度.依照索周法處理步驟,當(dāng)排列矩陣的列數(shù)近似等于數(shù)據(jù)所蘊(yùn)含的實(shí)際周期數(shù)值時(shí),各矩陣列數(shù)據(jù)的極值同時(shí)分布在不同矩陣的同一列中并使該列中存在較大的正數(shù)或較小的負(fù)數(shù),致使絕對(duì)值Δ較大.研究結(jié)果顯示:正、逆排列曲線波峰所對(duì)應(yīng)的矩陣列數(shù)即行星會(huì)合指數(shù)數(shù)值模擬數(shù)據(jù)所蘊(yùn)含的19.89 a顯著周期,并出現(xiàn)可能由蘊(yùn)含最顯著的19.89 a周期所疊合而成的新的39.35 a周期,同時(shí)也存在較明顯的12.93和25.20 a周期性規(guī)律.

圖3 1000—2000AD行星會(huì)合指數(shù)索周圖Fig.3 Cycle exploring method diagrams of planet juncture index from 1000 to 2000AD

分別繪制行星會(huì)合指數(shù)在1000—1500AD和1500—2000AD期間數(shù)值模擬數(shù)據(jù)的功率譜分析圖,如圖4所示.對(duì)圖譜進(jìn)行歸一化處理后不難發(fā)現(xiàn):兩段500 a等長(zhǎng)時(shí)間序列內(nèi),行星會(huì)合指數(shù)模擬數(shù)值均存在f=0.050 3 Hz的最強(qiáng)功率譜峰值,因此行星會(huì)合指數(shù)所凸顯的T=1/f=19.873 a周期最為顯著,這與采用索周法求得的行星會(huì)合指數(shù)方程周期性規(guī)律相吻合.劉復(fù)剛等[10]研究認(rèn)為,這一明顯的周期規(guī)律是由太陽(yáng)系中質(zhì)量最大的土木雙星運(yùn)動(dòng)軌道周期疊合造成的,推演的土木雙星會(huì)合周期也指向19.85 a,與行星會(huì)指數(shù)最明顯的19.873 a周期最為近似.對(duì)比圖3和圖4可以發(fā)現(xiàn):功率譜圖顯示的較強(qiáng)峰值所對(duì)應(yīng)的19.873、12.807 a周期與索周法處理結(jié)果中峰值所對(duì)應(yīng)的19.89、12.93 a周期存在顯著的一致性;同時(shí),索周法獲取的39.35 a周期可近似由兩倍的19.89 a顯著周期構(gòu)成.

圖4 1000—2000AD行星會(huì)合指數(shù)功率譜圖Fig.4 Power spectra of planet juncture index from 1000 to 2000AD

3 結(jié) 語(yǔ)

1) 行星會(huì)合指數(shù)運(yùn)動(dòng)學(xué)方程的創(chuàng)建,很大程度上將較為復(fù)雜的太陽(yáng)系內(nèi)行星之間的運(yùn)轉(zhuǎn)模式進(jìn)行簡(jiǎn)化,形象地展示出太陽(yáng)系內(nèi)行星系統(tǒng)會(huì)合和離散的動(dòng)態(tài)變化,推演出太陽(yáng)系行星系統(tǒng)復(fù)雜的繞轉(zhuǎn)運(yùn)動(dòng)周期.

2) 索周法在周期探求過(guò)程中凸顯出普適性、精確性等特點(diǎn),規(guī)避了典型周期處理方式涉及信號(hào)的連續(xù)性、基函數(shù)選擇、分解方式選擇的不確定性,較為全面準(zhǔn)確地分析出不同長(zhǎng)度數(shù)據(jù)所蘊(yùn)含的周期.根據(jù)周期分析原理和允許接受的誤差范圍,其結(jié)果精度由選取步長(zhǎng)決定,因此索周法所得結(jié)果在可控誤差范圍內(nèi)能很好地估算誤差值.考慮到實(shí)際計(jì)算的強(qiáng)度及可接受的探求結(jié)果精度,選取1 a為步長(zhǎng),所得結(jié)果的可控誤差范圍≤2.5%,完全適用于行星會(huì)合指數(shù)的周期性研究;若嘗試結(jié)合計(jì)算機(jī)參與索周分析過(guò)程,可將原先短周期性規(guī)律探求的優(yōu)勢(shì)進(jìn)一步覆蓋到中、長(zhǎng)周期規(guī)律探索中,進(jìn)一步彌補(bǔ)時(shí)間序列對(duì)于精準(zhǔn)長(zhǎng)周期求取的局限性.

3) 通過(guò)索周法對(duì)行星會(huì)合指數(shù)運(yùn)動(dòng)學(xué)方程在近千年時(shí)間尺度的周期性研究,分析出蘊(yùn)含的19.89 a顯著周期;同時(shí),存在可能由兩倍顯著的19.89 a周期構(gòu)成的39.35 a周期.這不僅為深入研究行星系統(tǒng)運(yùn)動(dòng)周期性規(guī)律開辟了新思路,也為研究行星周期性運(yùn)動(dòng)引發(fā)全球氣候及地表環(huán)境變化提供新的途徑和方法.