有限域上稀疏多元多項(xiàng)式插值算法*

唐 敏，鄧國(guó)強(qiáng)

桂林電子科技大學(xué) 數(shù)學(xué)與計(jì)算科學(xué)學(xué)院 廣西密碼學(xué)與信息安全重點(diǎn)實(shí)驗(yàn)室，廣西 桂林 541004

1 引言

插值是函數(shù)逼近的重要方法，利用它可通過(guò)函數(shù)在有限個(gè)點(diǎn)處的取值情況，估算出函數(shù)在其他點(diǎn)處的近似值。在信號(hào)處理[1]、不確定性量化[2-4]、壓縮感知[5]、圖像處理[6-7]等領(lǐng)域需要有效地利用較少的函數(shù)值恢復(fù)一個(gè)函數(shù)，而這個(gè)函數(shù)往往具有稀疏的表示形式，稀疏插值能有效地處理此類問(wèn)題。稀疏多元多項(xiàng)式插值是一種降低計(jì)算機(jī)代數(shù)算法時(shí)間復(fù)雜度的有效方法，是稀疏有理函數(shù)插值[8]及稀疏隱函數(shù)插值[9]的基本步驟，應(yīng)用在結(jié)式計(jì)算[9]、多元多項(xiàng)式GCD計(jì)算[10-12]、方程求根、數(shù)值微分、數(shù)值積分、曲面的外形設(shè)計(jì)和有限元法及上述提及的諸多應(yīng)用領(lǐng)域。因此提高稀疏多元多項(xiàng)式插值的效率十分有意義。

對(duì)于多元多項(xiàng)式，稠密插值所需的插值點(diǎn)個(gè)數(shù)為(d+1)n，其中d是變?cè)拇螖?shù)，n是變?cè)膫€(gè)數(shù)。研究工作表明，多元多項(xiàng)式插值算法的計(jì)算復(fù)雜度與變?cè)獋€(gè)數(shù)可以不呈指數(shù)關(guān)系，而是與目標(biāo)多項(xiàng)式的稀疏性相關(guān)，即是關(guān)于n、d、t和lbp的多項(xiàng)式函數(shù)，其中t是目標(biāo)多項(xiàng)式的項(xiàng)數(shù)，素?cái)?shù)p是有限域的特征。

1979年Zippel提出了第一個(gè)具有多項(xiàng)式時(shí)間復(fù)雜度的稀疏多元多項(xiàng)式插值算法[12]。該算法基于這樣的假設(shè)：如果在一個(gè)隨機(jī)的賦值點(diǎn)處多項(xiàng)式的取值為零，那么它就是一個(gè)零多項(xiàng)式。由于這個(gè)結(jié)論的成立是高概率的，因此Zippel算法是一個(gè)概率性的方法。在執(zhí)行上，Zippel算法順序地對(duì)每個(gè)變?cè)鹨徊逯?，不能并行化，插值點(diǎn)個(gè)數(shù)為ndt，時(shí)間復(fù)雜度為O(ndt3)。1990年，Zippel改進(jìn)了他的方法[13]，使用了形如的插值點(diǎn)，原算法中需要計(jì)算的一般線性系統(tǒng)變成了轉(zhuǎn)置的Vandermonde系統(tǒng)，求解時(shí)間從O(t3)降到O(t2)。由于Zippel算法是一種較好的具有多項(xiàng)式時(shí)間復(fù)雜度的方法，大部分計(jì)算機(jī)代數(shù)系統(tǒng)都內(nèi)置了基于Zippel算法的多元GCD計(jì)算的實(shí)現(xiàn)，比如Macsyma、Magma和Mathematica。

1988年，Ben-Or和Tiwari提出了一個(gè)確定性的算法[14]，可以對(duì)以整數(shù)、有理數(shù)、實(shí)數(shù)或復(fù)數(shù)為系數(shù)的多元多項(xiàng)式進(jìn)行稀疏插值。給定項(xiàng)數(shù)為t的多項(xiàng)式f的項(xiàng)數(shù)界T≥t，算法需要2t個(gè)插值點(diǎn)，即前n個(gè)素?cái)?shù) 的 冪 次，0≤i≤2t-1。該方法與Zippel算法不同，不是逐個(gè)變?cè)牟逯?，而是一次性地?duì)f的各個(gè)單項(xiàng)式進(jìn)行插值，而且可以并行化。主要的不足是插值點(diǎn)的存儲(chǔ)空間為tlbn位，當(dāng)多項(xiàng)式規(guī)模較大時(shí)，由于中間表達(dá)式膨脹而導(dǎo)致計(jì)算速度非常慢。

2000年，Kaltofen等設(shè)計(jì)了一個(gè)Zippel和Tiwari方法相結(jié)合的混合算法，稱之為“racing算法”[15-16]。為了減少插值點(diǎn)個(gè)數(shù)，在Zippel算法中插值下一個(gè)變?cè)臅r(shí)候，同時(shí)執(zhí)行Newton插值和單變?cè)腂en-Or/Tiwari插值，取其快者。

2010年，Javadi和Monagan提出了一個(gè)有限域上的稀疏多元多項(xiàng)式插值算法[17]，它也是一個(gè)概率算法。該方法是Ben-Or/Tiwari算法的一個(gè)改進(jìn)，使用O(t)個(gè)插值點(diǎn)獨(dú)立完成每個(gè)變?cè)牟逯?，共需O(nt)個(gè)插值點(diǎn)，并且可以并行化。為了準(zhǔn)確計(jì)算每個(gè)變?cè)诟鱾€(gè)單項(xiàng)式中的準(zhǔn)確次數(shù)，使用了二部圖完美匹配進(jìn)行測(cè)試和驗(yàn)證。該算法和Zippel算法及Kaltofen的racing算法進(jìn)行了比較，在運(yùn)行時(shí)間和插值點(diǎn)個(gè)數(shù)等方面具有一定的優(yōu)勢(shì)。

2011年，Cuyt和Lee給出了稀疏有理函數(shù)插值算法[8]，其中涉及到單變?cè)欣砗瘮?shù)插值和基于Zippel算法及Tiwari算法的稀疏多元多項(xiàng)式插值。2016年，基于Cuyt的稀疏有理函數(shù)插值算法，Tang給出了稀疏隱函數(shù)插值算法[9]。在稀疏有理函數(shù)插值和隱函數(shù)插值算法中，稀疏多元多項(xiàng)式插值作為基本的子函數(shù)，對(duì)整個(gè)算法的效率起到關(guān)鍵作用。

本文提出了一個(gè)有限域上的稀疏多元多項(xiàng)式插值算法，它是Javadi/Monagan算法的改進(jìn)，改進(jìn)方案主要有兩點(diǎn)：一是消除了必須給定項(xiàng)數(shù)界T的限制，通過(guò)計(jì)算預(yù)先設(shè)計(jì)好的特定矩陣的行列式，得到f的準(zhǔn)確項(xiàng)數(shù)；二是消除了必須給定次數(shù)界D的限制，通過(guò)構(gòu)造輔助函數(shù)，利用概率法結(jié)合提前終止技術(shù)的Cauchy插值法，計(jì)算出f的準(zhǔn)確次數(shù)，解決了Javadi和Monagan論文中提到的次數(shù)界D過(guò)高而導(dǎo)致高計(jì)算復(fù)雜度的問(wèn)題。更進(jìn)一步，改進(jìn)算法無(wú)需給定T和D，對(duì)于實(shí)際問(wèn)題在利用插值恢復(fù)或近似時(shí)更具實(shí)用性。

本文后續(xù)組織結(jié)構(gòu)如下：第2章介紹了Javadi/Monagan算法的基本思想；第3章給出了改進(jìn)的Javadi/Monagan算法；第4章分析了改進(jìn)算法的時(shí)間復(fù)雜度；第5章中設(shè)計(jì)了三組數(shù)值實(shí)驗(yàn)，對(duì)改進(jìn)算法和Javadi/Monagan算法進(jìn)行了比較，給出了實(shí)驗(yàn)結(jié)果；第6章給出了一個(gè)應(yīng)用實(shí)例。

2 Javadi/Monagan算法

本章給出稀疏多元多項(xiàng)式插值問(wèn)題的形式化定義及Javadi/Monagan算法的思想。

2.1 稀疏多元多項(xiàng)式插值問(wèn)題的形式化描述

令p是一個(gè)素?cái)?shù)，f∈Zp[x1，x2，…，xn]是一個(gè)用黑盒表示的多元多項(xiàng)式：

稀疏多元多項(xiàng)式插值的目標(biāo)是用盡可能少的插值點(diǎn)及較低的多項(xiàng)式時(shí)間復(fù)雜度算法還原多項(xiàng)式f。

2.2 Javadi/Monagan算法的思想

Javadi/Monagan算法可分成兩個(gè)步驟，首先確定各個(gè)單項(xiàng)式，然后確定系數(shù)。在第一個(gè)步驟中，逐一計(jì)算每個(gè)變?cè)诟鱾€(gè)單項(xiàng)式中的次數(shù)。第二個(gè)步驟需要求解一個(gè)Vandermonde線性代數(shù)系統(tǒng)。同時(shí)要求給定f的項(xiàng)數(shù)界T≥t及次數(shù)界D≥degf。

（1）確定單項(xiàng)式Mi，i=1，2，…，t。

③使用Berlekamp/Massey算法[18]生成n+1個(gè)關(guān)于z的單變?cè)囗?xiàng)式，使得，其中i=0，1，…，2t-1,j=0，1，…，n。

④確定變?cè)獂k，k=1，2，…，n在f中的單項(xiàng)式Mi，i=1，2，…，t中的次數(shù)eij。詳見2.3節(jié)。

⑤由④可確定eij，i=1，2，…，t，j=1，2，…，n，則單項(xiàng)式

（2）計(jì)算系數(shù)。

通過(guò)求解線性方程組計(jì)算系數(shù)ai，i=1，2，…，t。令r1，r2，…，rt是Λ1(z)的根，且ri=Mi(α1，α2，…，αn)，回顧是輸入為時(shí)黑盒的輸出，因此有：

此為Vandermonde系統(tǒng)，有唯一解。

2.3 確定變?cè)獂k的次數(shù)

令R1={r1，r2，…，rt}和分別是Λ1和Λk+1的全部根構(gòu)成的集合。令Dj包含單項(xiàng)式Mj中xk的所有可能的次數(shù)，即：

二部圖Gk定義如下：U和V是二部圖Gk中頂點(diǎn)個(gè)數(shù)為t的互補(bǔ)頂點(diǎn)子集，U和V中的結(jié)點(diǎn)表示R1和Rk中的元素，即ui∈U，vj∈V，分別用ri和標(biāo)記，ui和vj之間有一條權(quán)值為dij的邊當(dāng)且僅當(dāng)

引理1[17]變?cè)獂k在所有單項(xiàng)式中的次數(shù)可唯一確定當(dāng)且僅當(dāng)二部圖Gk有唯一的完美匹配。

以下給出為計(jì)算xk的次數(shù)，構(gòu)造的二部圖Gk沒有唯一的完美匹配時(shí)的解決方案。

隨機(jī)選擇不同于α1，α2，…，αn+1的元素令，使用 Berlekamp/Massey算法，以v′i=f(β′i)，0≤i≤ 2t-1為輸入，生成多項(xiàng)式Λ′k+1(z)。令是Λ′k+1(z)的根構(gòu)成的集合，G′k是通過(guò)Λ1(z)和Λ′k+1(z)構(gòu)造的二部圖。

定義1[17]二部圖G′k與Gk的交集定義如下：與G′k具有相同的結(jié)點(diǎn)，ri和之間有權(quán)值為dij的邊當(dāng)且僅當(dāng)在二部圖Gk中，ri連接且在二部圖G′k中，ri連接權(quán)值均為dij。

引理2[17]令eij=degxj(Mi)，在二部圖中，結(jié)點(diǎn)ri和之間有邊相連，且權(quán)值為eij。

定理1[17]令的交集有唯一的完美匹配。

引理1、引理2和定理1的證明參考文獻(xiàn)[17]。由定理1，如果Gk沒有唯一的完美匹配，可通過(guò)增加2t個(gè)插值點(diǎn)，構(gòu)造二部圖來(lái)確定變?cè)獂k在單項(xiàng)式Mj(1≤j≤t)中的次數(shù)。

3 改進(jìn)的Javadi/Monagan算法

本章給出對(duì)Javadi/Monagan算法進(jìn)行改進(jìn)的兩個(gè)主要策略。3.1節(jié)和3.2節(jié)分別描述了如何只利用黑盒的輸入輸出（即有限個(gè)插值點(diǎn)），計(jì)算多項(xiàng)式f的準(zhǔn)確項(xiàng)數(shù)t和準(zhǔn)確次數(shù)d，從而消除了原算法必須預(yù)先給定這兩個(gè)指標(biāo)的上界的限制，同時(shí)提高了實(shí)用性及降低了計(jì)算時(shí)間復(fù)雜度(分析過(guò)程見第4章)。

3.1 f的項(xiàng)數(shù)t的計(jì)算

除了Javadi/Monagan算法，許多插值算法都要求給定目標(biāo)多項(xiàng)式f的項(xiàng)數(shù)界T，本文給出一個(gè)準(zhǔn)確計(jì)算多項(xiàng)式f的項(xiàng)數(shù)t的方法。該方法由定理2保證。

定理2[14]令V是元素為(V)ij=vi+j-2的矩陣。Vl是由V的前l(fā)行前l(fā)列組成的方陣。如果t是多項(xiàng)式f的準(zhǔn)確項(xiàng)數(shù)，即f中非零單項(xiàng)式的個(gè)數(shù)，那么：

根據(jù)定理2，通過(guò)下述過(guò)程可計(jì)算f的準(zhǔn)確項(xiàng)數(shù)：

例如：黑盒多元多項(xiàng)式f定義為：令p=1 009，隨機(jī)生成中的3個(gè)整數(shù)，比如α1=66，α2=12，α3=3。

因此多項(xiàng)式f的準(zhǔn)確項(xiàng)數(shù)為5。

3.2 f的次數(shù)d的計(jì)算

在2.3節(jié)構(gòu)造二部圖的過(guò)程中，計(jì)算單項(xiàng)式Mj中xk的所有可能的次數(shù)Dj時(shí)，需要檢測(cè)的次數(shù)為D+1。因此Javadi/Monagan算法在給定的次數(shù)界過(guò)高時(shí)，會(huì)導(dǎo)致高計(jì)算復(fù)雜度。本節(jié)利用單變?cè)欣砗瘮?shù)插值計(jì)算f的準(zhǔn)確次數(shù)d，使得在確定每個(gè)變?cè)拇螖?shù)時(shí)，時(shí)間復(fù)雜度達(dá)到最低。方法描述如下：

（1）構(gòu)造輔助有理函數(shù)H。

通過(guò)在f中引入齊次變?cè)獄，構(gòu)造f為分子，g為分母的輔助有理函數(shù)H。其中f(x1，x2，…，xn)為黑盒多項(xiàng)式，g為隨機(jī)生成的關(guān)于變?cè)獄的一次多項(xiàng)式。若視F為關(guān)于變?cè)獄的單變?cè)囗?xiàng)式，則系數(shù)為關(guān)于x1，x2，…，xn的多元多項(xiàng)式，注意到系數(shù)多項(xiàng)式與z是齊次的。

易證f的全次數(shù)與F中z的最高次相同。

（2）對(duì)H進(jìn)行單變?cè)欣砗瘮?shù)插值。

任取n元組 (α1，α2，…，αn)，對(duì)H進(jìn)行關(guān)于變?cè)獄的單變?cè)欣砗瘮?shù)插值，可得到：

（3）H(α1，α2，…，αn，z)的分子f(α1z，α2z，…，αnz)關(guān)于z的最高次即為黑盒多項(xiàng)式f的全次數(shù)。

假設(shè)f(α1z，α2z，…，αnz)和g(z)有非平凡的公因子，即g(z)|f(α1z，α2z，…，αnz)，也就是說(shuō)，g(z)的根也是f(α1z，α2z，…，αnz)的根。因?yàn)間(z)是隨機(jī)生成的，而f(α1z，α2z，…，αnz)的根是有窮的，所以f(α1z，α2z，…，αnz)和g(z)高概率互素，從而f(α1z，α2z，…，αnz)關(guān)于z的次數(shù)即為黑盒多項(xiàng)式f的全次數(shù)。

令 (α1，α2)=(1，1)，g=5z+3 ，則：

利用單變?cè)欣砗瘮?shù)插值可得到H(1，1，z)=(z4+2z2)/(5z+3)，則f的準(zhǔn)確次數(shù)為H的分子中z的最高次，即為4。

接下來(lái)給出單變?cè)欣砗瘮?shù)插值的實(shí)現(xiàn)過(guò)程。對(duì)于單變?cè)欣砗瘮?shù)f/g∈Κ(Z)，應(yīng)用文獻(xiàn)[19]提出的概率方法，并結(jié)合提前終止技術(shù)的Cauchy插值來(lái)恢復(fù)f/g。來(lái)源于文獻(xiàn)[19-20]的引理給出了單變?cè)欣砗瘮?shù)插值的理論基礎(chǔ)。

引理3[20]令Κ是任意域，F(xiàn)(Z)，G(Z)，H(Z)∈Κ(Z)且 gcd(F，G)=1。令是非負(fù)整數(shù)，且1。令ik，f/g∈Κ(Z)不必是Κ中互不相同的元素，滿足：

對(duì)于l≥ 2 ，令hl(Z)，ql(Z)∈Κ(Z)分別是Euclidean多項(xiàng)式余項(xiàng)序列中的第l次余項(xiàng)和商。

對(duì)于l≥ 2 ，令wl(Z)，gl(Z)∈Κ(Z)是擴(kuò)展Euclidean方法中的乘子wlh0+glh1=hl，即：

那么存在下標(biāo)j≥ 1滿足，并且對(duì)下標(biāo)j有：

根據(jù)引理3，如果給定單變?cè)欣砗瘮?shù)的次數(shù)上界，使用擴(kuò)展Euclidean算法可計(jì)算有理函數(shù)的分子和分母。在文獻(xiàn)[19]中，提出了基于引理3的單變?cè)欣砗瘮?shù)插值算法。為確定準(zhǔn)確的f/g的次數(shù)，可對(duì)k從1迭代到d+e+1，其中d、e分別是f和g的次數(shù)。

基于概率法和提前終止技術(shù)的Cauchy插值法恢復(fù)單變?cè)欣砗瘮?shù)，算法描述如下：

算法1單變?cè)欣砗瘮?shù)插值算法

其 中 Interpolate(i1，i2，…，ik;H(i1)，H(i2)，…，H(ik)) 表示給定k個(gè)點(diǎn) (i1，H(i1))，(i2，H(i2))，…，(ik，H(ik))的自變量和函數(shù)值，利用Newton插值或Lagrange插值得到次數(shù)為k-1的單變?cè)囗?xiàng)式；Rem表示余式，Quo表示商。

黑盒多項(xiàng)式f的次數(shù)即為deg(h1)。

3.3 改進(jìn)算法

本節(jié)給出改進(jìn)的Javadi/Monagan算法的偽代碼。

算法2改進(jìn)的Javadi/Monagan算法(有限域上稀疏多元多項(xiàng)式插值算法)

以下給出改進(jìn)的稀疏多元多項(xiàng)式插值算法的一個(gè)實(shí)例。黑盒多元多項(xiàng)式f定義為：

令p=101，變?cè)獋€(gè)數(shù)n=3。重構(gòu)f的過(guò)程如下：

（2）計(jì)算出多項(xiàng)式f的準(zhǔn)確項(xiàng)數(shù)t=5。

（4）計(jì)算Λk(z)，k=1，2，…，4:

令vi=B(βi)，i=0，1，…，2t-1，當(dāng)k＞1時(shí)，用替換使用Berlekamp/Massey算法生成由序列v0，v1，…，v2t-1得到的Λk∈Zp(z)：

（5）計(jì)算出多項(xiàng)式f的次數(shù)d=5。

（6）確定degxk(Mi),1≤i≤t,k=1，2，…，n。

以第1個(gè)變?cè)獂1為例，計(jì)算Λ1的根：

計(jì)算Λ2的根：

測(cè)試rj×(αn+1/α1)i，0≤i≤d，j=1，2，…，t,得到：

構(gòu)造二部圖G1，如圖1所示。

Fig.1 Bipartite graphG1圖1 二部圖G1

因?yàn)槎繄DG1沒有唯一的完美匹配，所以生成α5=54。

計(jì)算Λ5(z)=z5+100z4+73z3+85z2+51z+94及Λ5的根。構(gòu)造二部圖G1′，如圖2所示。

Fig.2 Bipartite graphG1′圖2 二部圖G1′

構(gòu)造Gk和Gk′的交集，如圖3所示。因此，x1在Mi,1≤i≤t中的次數(shù)分別為0,0,2,2,0。同理，計(jì)算x2和x3在Mi中的次數(shù)，分別為0,0,1,2,1和0,5,1,1,2。則

Fig.3 Bipartite graph圖3 二部圖

4 時(shí)間復(fù)雜度分析

本章討論有限域上改進(jìn)的稀疏多元多項(xiàng)式插值算法（算法2）的時(shí)間復(fù)雜度。需要考慮的主要因素有：

（1）計(jì)算多項(xiàng)式f的項(xiàng)數(shù)t，時(shí)間復(fù)雜度為O(t2)。

（2）計(jì)算多項(xiàng)式f的次數(shù)d，時(shí)間復(fù)雜度為O(d)。

（3）調(diào)用Berlekamp/Massey算法計(jì)算Λk至多n+2次，每次調(diào)用的復(fù)雜度為O(t2)。

（4）使用文獻(xiàn)[13]中的技術(shù)求解Vandermonde系統(tǒng)，時(shí)間復(fù)雜度為O(t2)。

（5）求解Λk的根，時(shí)間復(fù)雜度為O(t2lbp)。

（6）計(jì)算變?cè)獂k的次數(shù)，時(shí)間復(fù)雜度為O(tlbt+dtlbt)。

（7）構(gòu)造二部圖Gk，時(shí)間復(fù)雜度為O(dt2)，求Gk和G′k的交集，時(shí)間復(fù)雜度為O(dtlbd)。

綜上所述，改進(jìn)算法的時(shí)間復(fù)雜度為：

從時(shí)間復(fù)雜度分析可以看出，如果給定的項(xiàng)數(shù)界T和次數(shù)界D過(guò)高，會(huì)導(dǎo)致高計(jì)算復(fù)雜度，而改進(jìn)算法花費(fèi)較少的代價(jià)計(jì)算出多項(xiàng)式f的準(zhǔn)確項(xiàng)數(shù)t及次數(shù)d，達(dá)到了此類型算法在后續(xù)操作上的最低時(shí)間復(fù)雜度。

5 數(shù)值實(shí)驗(yàn)

本章對(duì)改進(jìn)算法和Javadi/Monagan算法進(jìn)行性能比較。兩種算法的編程環(huán)境均為Maple15，程序運(yùn)行的硬件環(huán)境為Intel?CoreTMi7 2.20 GHz處理器和4.00 GB內(nèi)存，操作系統(tǒng)為Windows 7。注意兩個(gè)程序都是順序執(zhí)行，在確定變?cè)獂1，x2，…，xn的次數(shù)時(shí)未并行化。

本章給出了3組測(cè)試集的性能比較結(jié)果，使用的多項(xiàng)式都是隨機(jī)生成的，比較的對(duì)象為CPU時(shí)間。黑盒中的多元多項(xiàng)式系數(shù)取自Zp，其中p=100 003。

測(cè)試集1本組測(cè)試集為6個(gè)包含3個(gè)變?cè)亩嘣囗?xiàng)式。第i個(gè)多項(xiàng)式（1≤i≤6）使用如下的Maple命令隨機(jī)生成：

其中，第i個(gè)多項(xiàng)式包含2i個(gè)非零項(xiàng)，d=20是多項(xiàng)式的準(zhǔn)確次數(shù)。該組測(cè)試分別執(zhí)行了改進(jìn)算法和Javadi/Monagan算法，記錄每個(gè)多項(xiàng)式在兩種算法下的運(yùn)行時(shí)間（單位：s），如表1所示。表頭第1列Ex.表示例子的編號(hào)；第2列t表示項(xiàng)數(shù)；第3列括號(hào)中的百分比表示改進(jìn)算法計(jì)算項(xiàng)數(shù)和次數(shù)的總時(shí)間與整個(gè)算法的運(yùn)行時(shí)間的比值；改進(jìn)算法無(wú)需給定項(xiàng)數(shù)界和次數(shù)界，Javadi/Monagan算法執(zhí)行時(shí)，分別給定的次數(shù)界為D=20，D=30和D=30，項(xiàng)數(shù)界為多項(xiàng)式的準(zhǔn)確項(xiàng)數(shù)T=2i(1≤i≤6)。表2和表3的設(shè)置同表1。

Table 1 Performance comparison of improved algorithm and Javadi/Monagan algorithm(n=3)表1 改進(jìn)算法與Javadi/Monagan算法的比較(n=3)

從表1可以看出，隨著i的增加，兩種算法的執(zhí)行時(shí)間也隨之增加。在給定準(zhǔn)確項(xiàng)數(shù)界T=2i和準(zhǔn)確次數(shù)界D=20的情況下，Javadi/Monagan算法的執(zhí)行時(shí)間優(yōu)于改進(jìn)算法，原因是改進(jìn)算法利用3.1節(jié)和3.2節(jié)的方法分別計(jì)算出準(zhǔn)確的項(xiàng)數(shù)和次數(shù)，不需預(yù)先給定這兩個(gè)數(shù)據(jù)。顯然，這兩者的計(jì)算需要耗費(fèi)一定的時(shí)間，但從表1中的第3列括號(hào)中的百分比可以看出，這兩者的計(jì)算時(shí)間占整體運(yùn)行時(shí)間的比值隨著次數(shù)的增加越來(lái)越小，因而優(yōu)勢(shì)越來(lái)越明顯。另一方面，隨著次數(shù)界D的增加，Javadi/Monagan算法需要的執(zhí)行時(shí)間也隨之增加，當(dāng)D較大時(shí)，算法的時(shí)間復(fù)雜度較高，這就是Javadi和Monagan在論文中提到的壞次數(shù)界問(wèn)題（bad degree bound）。改進(jìn)算法則有效地避免了這一問(wèn)題。

測(cè)試集2本組測(cè)試集為6個(gè)包含6個(gè)變?cè)亩嘣囗?xiàng)式。第i個(gè)多項(xiàng)式(1≤i≤6)使用如下的Maple命令隨機(jī)生成：

表2給出了測(cè)試集2下改進(jìn)算法和Javadi/Monagan算法的執(zhí)行時(shí)間，表格的各項(xiàng)含義與測(cè)試集1相同。

測(cè)試集3本組測(cè)試集為6個(gè)包含12個(gè)變?cè)亩嘣囗?xiàng)式。第i個(gè)多項(xiàng)式(1≤i≤6)使用如下的Maple命令隨機(jī)生成：

Table 2 Performance comparison of improved algorithm and Javadi/Monagan algorithm(n=6)表2 改進(jìn)算法與Javadi/Monagan算法的比較(n=6)

表3給出了測(cè)試集3下改進(jìn)算法和Javadi/Monagan算法的執(zhí)行時(shí)間，表格的各項(xiàng)含義與測(cè)試集1相同。

Table 3 Performance comparison of improved algorithm and Javadi/Monagan algorithm(n=12)表3 改進(jìn)算法與Javadi/Monagan算法的比較(n=12)

從3組測(cè)試集可以看出，改進(jìn)算法在確定插值多項(xiàng)式f的項(xiàng)數(shù)和次數(shù)上，花費(fèi)的時(shí)間占整個(gè)算法運(yùn)行時(shí)間的比例非常小，多項(xiàng)式規(guī)模越大，比例越小。如果Javadi/Monagan算法執(zhí)行時(shí)，給定的次數(shù)界D越高，運(yùn)行時(shí)間越長(zhǎng)，而改進(jìn)算法不受此影響。

6 應(yīng)用實(shí)例

本章給出改進(jìn)的Javadi/Monagan稀疏插值算法在幾何問(wèn)題上的一個(gè)應(yīng)用。著名的Morley三等分定理描述如下：在任意三角形ABC中，對(duì)3個(gè)內(nèi)角∠CAB、∠ABC和∠ACB分別進(jìn)行三等分，相鄰邊交于三點(diǎn)P、Q、R，則三角形PQR必為一個(gè)正三角形，如圖4所示。令a=BC,b=CA,c=AB,x=PQ=QR=RP,求x和a、b、c之間的關(guān)系。

Fig.4 Morley triangle圖4 Morley三角形

根據(jù)題設(shè)中的幾何關(guān)系和相關(guān)定理可知（推導(dǎo)過(guò)程略）：

其中，R是三角形ABC的外接圓半徑，有如下方程組：

令x和a、b、c之間的關(guān)系為f(x，a，b，c)=0 ，如果利用結(jié)式等符號(hào)消元法，消去方程組（1）中的變?cè)猵、q、r，那么可得到x和a、b、c之間的關(guān)系f(x，a，b，c)。實(shí)際上，如果把a(bǔ)、b、c視為符號(hào)參數(shù)，消元過(guò)程因中間表達(dá)式膨脹而無(wú)法完成，若將a、b、c用實(shí)數(shù)α1、α2、α3進(jìn)行替換，消元過(guò)程可順利進(jìn)行，結(jié)果是實(shí)例化的關(guān)系式f(x，α1，α2，α3)，例如（假定模p=105）：

如果視f(x，α1，α2，α3)為關(guān)于變?cè)獂的單變?cè)囗?xiàng)式，則x的各次冪的系數(shù)是關(guān)于變?cè)猘、b、c的多元多項(xiàng)式，使用稀疏多元多項(xiàng)式插值恢復(fù)系數(shù)多項(xiàng)式，即可得到x和a、b、c之間的關(guān)系。在此例中，x的各個(gè)冪次的系數(shù)多項(xiàng)式的項(xiàng)數(shù)和次數(shù)都無(wú)法估計(jì)，如果采用原始的Javadi/Monagan算法，可能因?yàn)轫?xiàng)數(shù)界和次數(shù)界設(shè)置不正確而導(dǎo)致無(wú)法得出結(jié)果或計(jì)算復(fù)雜度過(guò)高。若采用改進(jìn)的Javadi/Monagan算法，可準(zhǔn)確計(jì)算出各個(gè)系數(shù)多項(xiàng)式的項(xiàng)數(shù)和次數(shù)，最終得到的x和a、b、c之間的關(guān)系式：

此例中，f(x，a，b，c)共有804項(xiàng)，其中x27的系數(shù)多項(xiàng)式的項(xiàng)數(shù)為45，次數(shù)為32。x25的系數(shù)多項(xiàng)式的項(xiàng)數(shù)為28，次數(shù)為34。其余系數(shù)多項(xiàng)式的指標(biāo)不再贅述。

在結(jié)式消元、GCD計(jì)算、幾何組合優(yōu)化、信號(hào)處理等問(wèn)題上，很多情況下無(wú)法預(yù)知目標(biāo)多項(xiàng)式的次數(shù)和項(xiàng)數(shù)，故改進(jìn)算法針對(duì)這些問(wèn)題更具實(shí)用性。

7 結(jié)束語(yǔ)

稀疏多元多項(xiàng)式插值是很多計(jì)算機(jī)代數(shù)算法的子函數(shù)，也廣泛應(yīng)用于信號(hào)處理、壓縮感知、圖像處理等領(lǐng)域。在給定項(xiàng)數(shù)界T和次數(shù)界D較高的情況下，改進(jìn)算法可有效降低Javadi/Monagan算法的計(jì)算復(fù)雜度。進(jìn)一步，改進(jìn)算法消除了傳統(tǒng)插值算法中必須預(yù)先給定插值多項(xiàng)式項(xiàng)數(shù)T和次數(shù)D的必要條件，更具實(shí)用性。

3.3節(jié)的算法2給出了有限域上改進(jìn)的Javadi/Monagan稀疏多元多項(xiàng)式插值算法，其中兩個(gè)循環(huán)體（算法2中步驟4到步驟8及步驟10到步驟19）的計(jì)算任務(wù)可并行化處理，一是生成關(guān)于變?cè)獄的單變?cè)囗?xiàng)式Λ1(z)，Λ2(z)，…，Λn+1(z)，可分為相互獨(dú)立的n+1個(gè)子任務(wù)，二是確定變?cè)獂1，x2，…，xn在各個(gè)單項(xiàng)式中的次數(shù)，可分為相互獨(dú)立的n個(gè)子任務(wù)。假設(shè)有q個(gè)計(jì)算內(nèi)核的集群系統(tǒng)，每個(gè)內(nèi)核可獨(dú)立承擔(dān)上述計(jì)算Λj(z)，j=1，2，…，n+1 的 (n+1)/q個(gè)子任務(wù)，以及計(jì)算變?cè)獂j，j=1，2，…，n在各個(gè)單項(xiàng)式Mi，i=1，2，…，t中的次數(shù)eij的n q個(gè)子任務(wù)，計(jì)算結(jié)果可保存在一個(gè)文件中，最終可確定黑盒多項(xiàng)式f中的各個(gè)單項(xiàng)式。在后續(xù)工作中，擬將改進(jìn)算法進(jìn)行并行化處理，以進(jìn)一步提高稀疏多元多項(xiàng)式插值算法的效率。