改進(jìn)的超越函數(shù)分段線性逼近方法

田征 杜慧敏 黃小康

摘要：針對(duì)超越函數(shù)計(jì)算中所采用的分段線性逼近算法存在的無(wú)法提前確定精度及部分區(qū)間資源浪費(fèi)的問(wèn)題，提出一種改進(jìn)的分段線性逼近超越函數(shù)算法。該算法由預(yù)定義的逼近區(qū)間端點(diǎn)計(jì)算出用于逼近的線性函數(shù)，根據(jù)被逼近函數(shù)的凹凸性對(duì)所計(jì)算線性函數(shù)進(jìn)行調(diào)整，在此基礎(chǔ)上計(jì)算出預(yù)定義逼近區(qū)間內(nèi)調(diào)整后函數(shù)與被逼近函數(shù)之間的最大誤差；按照所需精度的要求，自動(dòng)調(diào)整逼近區(qū)間，通過(guò)該過(guò)程的迭代，獲得了較少分段次數(shù)。算法結(jié)果在Matlab上進(jìn)行仿真，仿真結(jié)果表明，所提算法的分段數(shù)相比等分法減少了60%。所提算法在保證精度的前提下，降低了查找表（LUT）的資源消耗。

關(guān)鍵詞：

分段線性逼近；超越函數(shù)；查找表；資源浪費(fèi)；優(yōu)化分段方法

中圖分類號(hào)： TP391.75 文獻(xiàn)標(biāo)志碼：A

0引言

圖形處理器（Graphics Processing Unit， GPU）是各種嵌入式系統(tǒng)、個(gè)人機(jī)（Personal Computer， PC）、工作站和游戲機(jī)等不可缺少的重要部件。浮點(diǎn)超越函數(shù)單元是GPU數(shù)據(jù)通路中的重要部件，其性能直接影響圖形渲染效果[1-2]。

目前在現(xiàn)場(chǎng)可編程門陣列（FieldProgrammable Gate Array， FPGA）中計(jì)算超越函數(shù)常用的方法有級(jí)數(shù)近似法、查表法（Look Up Table， LUT）、坐標(biāo)旋轉(zhuǎn)數(shù)字計(jì)算（Coordinate Rotation Digital Computer， CORDIC）算法和分段線性逼近法等。

其中，級(jí)數(shù)近似法展開式較為復(fù)雜，硬件實(shí)現(xiàn)復(fù)雜，資源消耗較大。查表法雖然計(jì)算簡(jiǎn)單、易于實(shí)現(xiàn)，但是所需存儲(chǔ)單元隨著計(jì)算精度的提高呈指數(shù)形式增加，資源消耗大[3]。CORDIC算法[4]作為一種便于FPGA實(shí)現(xiàn)的超越函數(shù)計(jì)算方法得到了廣泛關(guān)注，但其收斂速度與數(shù)據(jù)的表示精度成反比，當(dāng)精度要求較高時(shí)，算法的迭代次數(shù)較多，計(jì)算延遲會(huì)增大。與之相比，分段線性逼近法[5-6]將低階多項(xiàng)式與較小的查找表相結(jié)合，資源消耗少，速度較快，被廣泛應(yīng)用于傳感網(wǎng)絡(luò)中的數(shù)據(jù)壓縮[7]、非線性模型到線性模型的轉(zhuǎn)換[8]以及圖形圖像處理[9]等領(lǐng)域，已成為在計(jì)算資源有限條件下超越函數(shù)計(jì)算的一種較為理想的選擇。

為了提高函數(shù)的計(jì)算速度，在硬件設(shè)計(jì)中一般利用查找表實(shí)現(xiàn)分段線性逼近算法[6]，段數(shù)越多，查找表越大，誤差則越小，因此，分段線性逼近算法需要在查找表、精度和分段數(shù)目之間尋求一個(gè)合理的平衡。其中有均方誤差法[10]、區(qū)間2k等分法[11]、面積法[12]等多種分段算法，但這些算法都是在計(jì)算完成后才可以知道分段精度，難以在精度已知的條件下完成分段數(shù)的計(jì)算。

文獻(xiàn)[6]中提出了一種超越函數(shù)計(jì)算的最佳等距分段線性逼近（Optimal EquiDistant PieceWise Linear approximation，OED_PWL）方法，該方法可以通過(guò)調(diào)節(jié)分段數(shù)目來(lái)完成計(jì)算精度的靈活控制，從而可以在保證精度的條件下完成分段數(shù)的計(jì)算。相較區(qū)間2k等分法其算法性能、資源消耗等方面的表現(xiàn)更好；但是，根據(jù)文獻(xiàn)[6]提出的方法進(jìn)行tanh函數(shù)計(jì)算時(shí)，部分區(qū)間的計(jì)算誤差較大，原因是tanh函數(shù)斜率在某些分段較平緩，而在另外一些分段斜率較陡峭，若采用均勻分段的方法，則存在分段平緩處的資源浪費(fèi)問(wèn)題。在后面的實(shí)驗(yàn)結(jié)果中，本文會(huì)與OED_PWL算法進(jìn)行對(duì)比。

因此本文提出一種改進(jìn)的超越函數(shù)分段線性逼近方法，通過(guò)自動(dòng)識(shí)別逼近誤差進(jìn)行分段使其不斷細(xì)分，實(shí)現(xiàn)了一個(gè)精度已知條件下較優(yōu)的分段方法。

1分段線性逼近算法基本原理

分段線性逼近算法的原理是把非線性特性曲線分成若干個(gè)區(qū)段，在每個(gè)區(qū)段中用直線段近似地逼近特性曲線。算法原理如圖1所示。

分段線性逼近算法是根據(jù)一定的分段方法，如面積法、等分法等，將曲線f（x）所在的一段區(qū)間[A，B]劃分為若干段，在每一段里面通過(guò)一個(gè)線性函數(shù)來(lái)進(jìn)行逼近：

y（x）ax+b（1）

其中：a表示線性函數(shù)的斜率；b表示線性函數(shù)的偏移量。其中，在每一段逼近區(qū)間[x1，x2]里面，區(qū)間的兩個(gè)端點(diǎn)需要滿足：

{f（x1）=y（x1）， f（x2）=y（x2）}（2）

在通過(guò)某種方法得到[A，B]上超越函數(shù)曲線的分段情況以及各個(gè)段的逼近系數(shù)ai、bi后，即可完成超越函數(shù)的分段線性逼近。

算法硬件實(shí)現(xiàn)分為3步：

步驟1按照某種方法進(jìn)行分段區(qū)間的劃分，求取所有分段點(diǎn)。

步驟2在每個(gè)分段區(qū)間，根據(jù)兩個(gè)端點(diǎn)的值求取直線段f（x）=ax+b，將每一段的ai和bi存放在查找表中以待使用。

步驟3利用查找表，讀取每一段直線方程的系數(shù)，通過(guò)乘法和加法計(jì)算每一段函數(shù)的逼近值，從而實(shí)現(xiàn)超越函數(shù)的分段線性逼近。

2改進(jìn)的分段線性逼近算法

2.1算法原理

如何在給定的精度下，得到一個(gè)較好的分段，使得應(yīng)用分段計(jì)算值滿足給定的精度要求是分段線性逼近算法的一個(gè)難題。本文提出一種基于區(qū)間2k等分法，通過(guò)迭代來(lái)獲得分段區(qū)間的分段逼近方法，接下來(lái)對(duì)其原理進(jìn)行詳細(xì)介紹。

區(qū)間2k等分法會(huì)將[A，B]上的曲線均勻地分為2k段，k的取值取決于查找表的大小及資源的分配情況。這是分段線性逼近中最常用的一種分段方法，具有實(shí)現(xiàn)簡(jiǎn)單、精度較高等優(yōu)點(diǎn)；但是由于超越函數(shù)曲線的斜率在定義域內(nèi)是不均勻分布的，部分計(jì)算域內(nèi)曲線段斜率較大，因此采用直線段線性逼近誤差較大；而另外一部分計(jì)算域曲線段斜率較小，逼近誤差較小。為了滿足誤差較大曲線段的誤差要求，往往會(huì)導(dǎo)致分段數(shù)過(guò)多，導(dǎo)致斜率較小曲線段的資源浪費(fèi)。

因此本文提出一種改進(jìn)型的超越函數(shù)分段線性逼近算法，可以在精度已知的前提下，完成對(duì)分段數(shù)的最優(yōu)選取，使得計(jì)算函數(shù)誤差在一定的范圍內(nèi)。算法的主要思想為：對(duì)斜率較大、誤差較大的函數(shù)區(qū)間進(jìn)行進(jìn)一步細(xì)分，避免了對(duì)曲線平滑部分分段的浪費(fèi)，減少了查找表資源的消耗。

以對(duì)數(shù)lb x在區(qū)間[1，5]的計(jì)算說(shuō)明算法的思想：首先整個(gè)逼近區(qū)間等分為[1，3]和[3，5]回復(fù)：分段區(qū)間沒(méi)有問(wèn)題兩段，并對(duì)兩段分別進(jìn)行線性逼近，得到line1的兩條直線段，可以看出，在[1，3]區(qū)間的直線逼近誤差比較大，因此，對(duì)區(qū)間[1，3]繼續(xù)進(jìn)行細(xì)分，分為更小的兩個(gè)區(qū)間：[1，2]和[2，3]在等于2時(shí)，兩個(gè)區(qū)間都包含在內(nèi)，其中一個(gè)需用開區(qū)間來(lái)表示，請(qǐng)明確?；貜?fù)：分段區(qū)間沒(méi)有問(wèn)題，并對(duì)其分別進(jìn)行線性逼近，得到兩條誤差較低的直線段line2，完成分段區(qū)間的細(xì)分過(guò)程。

2.2算法實(shí)現(xiàn)流程

本文實(shí)現(xiàn)的算法流程如圖3所示，通過(guò)迭代的方法來(lái)完成分段區(qū)間的計(jì)算。

圖4說(shuō)明了算法的一次逼近過(guò)程，其中逼近區(qū)間段的中點(diǎn)處逼近誤差為d，假設(shè)原逼近直線y0=ax+b向曲線方向偏移d/2，得到新的逼近直線y1=ax+b+d/2，新的直線更加逼近曲線。

本文在求取被逼近函數(shù)f（x）與新的逼近直線y1 （x）之間的最大誤差dmax=f（x）-y1（x），對(duì)dmax求導(dǎo)數(shù)，令導(dǎo)數(shù)等于零，將得到x的值，代入dmax計(jì)算出最大誤差。

在本文的改進(jìn)型分段線性逼近算法中，所有數(shù)據(jù)采用雙精度類型，在每個(gè)分段區(qū)間上，判斷被逼近函數(shù)的凹凸性以及該函數(shù)在區(qū)間內(nèi)是否存在拐點(diǎn)，求取逼近直線以及被逼近函數(shù)之間的最大誤差，然后與給定誤差進(jìn)行比較：若小于給定誤差，則停止該段的繼續(xù)細(xì)分；若大于給定誤差，則對(duì)該段逼近區(qū)間繼續(xù)細(xì)分。

若刪除，可從此處刪除偽代碼。

下面是改進(jìn)型的分段線性逼近算法的C語(yǔ)言偽代碼。

圖3的表述，與此處的“c語(yǔ)言偽代碼”？是否有重復(fù)表達(dá)現(xiàn)象，是否可以刪除某一種表述方式？請(qǐng)明確。

回復(fù)：其中分段區(qū)間沒(méi)有問(wèn)題，還有就是那個(gè)關(guān)于偽代碼和流程圖重復(fù)問(wèn)題，我和我們老師商量了，可以刪除也可以不刪除，如果要?jiǎng)h除，就把偽代碼部分刪除。（附件是刪除偽代碼的稿件）

：

程序前

void subsection（double x0， double x1）

{

double y0=x0坐標(biāo)處要逼近的函數(shù)值f（x0）；

double y1=x1坐標(biāo)處要逼近的函數(shù)值f（x1）；

y′=f ″（x）=sign/*對(duì)逼近函數(shù)求導(dǎo)判斷凹凸性，sign=1表示凸函數(shù)，sign=0表示凹函數(shù)*/

令f ′（c）=0求出c值，c則為函數(shù)逼近區(qū)間內(nèi)的拐點(diǎn)；

if（y′=0）/*函數(shù)有拐點(diǎn)*/

{

重新劃分逼近區(qū)間為[x0，c]及[c，x1]，確保每段新的逼近區(qū)間都是單調(diào)的；

}

else；

mid=（x0+x1）/2.0；

a=（y1-y0）/（x1-x0）；

b=y1-a*（x1）；

max=d=f（mid）-（a*mid+b）；/*mid坐標(biāo)處的逼近誤差*/

if（sign==1）

b1=b+d/2.0；

else

b1=b-d/2.0；

y1（x）=ax+b1；//得到新的逼近直線y1

d=f（x）-y//被逼近函數(shù)與逼近函數(shù)的差值

d′=f ′（x）-y′//對(duì)函數(shù)d求導(dǎo)數(shù)

d′=f（m）-y1′（m）=0//求函數(shù)的極值點(diǎn)m

max=f（m）-y1（m）；//逼近區(qū)間內(nèi)的最大誤差

if （max<指定誤差）{

打印x0，x1坐標(biāo)，并打印逼近直線參數(shù)m及n到指定文檔；

}

else {

subsection（x0，mid）；

subsection（mid，x1）；

}

}

//說(shuō)明：x0與x1為要逼近的函數(shù)取值范圍；

程序后

3實(shí)驗(yàn)結(jié)果與對(duì)比

為了更好地說(shuō)明問(wèn)題，采用1E-3和1E-4作為預(yù)設(shè)的精度，可以滿足移動(dòng)端圖形圖像處理等應(yīng)用的精度要求；取值范圍設(shè)定為[1，2]，可以充分地說(shuō)明逼近的正確性。表1為各個(gè)超越函數(shù)在區(qū)間2k等分法和本文算法下得到的分段數(shù)量對(duì)比。

可以看出，改進(jìn)后的分段線性逼近算法對(duì)斜率變化較大的超越函數(shù)作用會(huì)比較明顯。而且如果精度提升，2k分段數(shù)會(huì)劇烈增長(zhǎng)，而本文提出的算法分段數(shù)增加較慢。

表2為本文算法與其他兩種算法在相同實(shí)驗(yàn)條件下的分段結(jié)果比較，可以看出：在保證誤差精度的前提下，本文算法的分段數(shù)相比等分法減少了60%以上，相比OED_PWL法減少了20%。對(duì)于FPGA硬件實(shí)現(xiàn)而言，分段數(shù)越少，需要的查找表則越小，從而節(jié)省了硬件查找表資源。

針對(duì)超越函數(shù)中的對(duì)數(shù)函數(shù)lb （1+x），文獻(xiàn)[9，13-14]等多篇文獻(xiàn)都對(duì)其進(jìn)行了分段區(qū)間的優(yōu)化。表3為本文算法對(duì)對(duì)數(shù)函數(shù)的分段逼近結(jié)果與其他文獻(xiàn)算法的對(duì)比，由于其他文獻(xiàn)的分段算法沒(méi)有做到精度已知下的分段劃分，因此實(shí)現(xiàn)精度各不相同。

由表2與表3的實(shí)驗(yàn)結(jié)果可以得出：

1）與其他已經(jīng)提出的算法相比，在相同精度下，本文的算法劃分的逼近區(qū)間更少，則查找到所屬分段的時(shí)間就會(huì)變少，所占據(jù)的存儲(chǔ)空間更小，更適合用硬件實(shí)現(xiàn)。

2）根據(jù)表3，與文獻(xiàn)[9]相比，分段數(shù)目相差不多的情況下，本文的算法可以大大降低逼近誤差。

綜上所述，改進(jìn)型的分段線性逼近算法同樣適用于tanh（x）函數(shù)、lb x函數(shù)等多種超越函數(shù)算法。

4對(duì)數(shù)函數(shù)的FPGA實(shí)現(xiàn)

完成算法的基本設(shè)計(jì)之后，本文將其應(yīng)用在某個(gè)超越函數(shù)運(yùn)算單元中，用于計(jì)算以2為底的對(duì)數(shù)轉(zhuǎn)換結(jié)果，其電路結(jié)構(gòu)如圖5所示。

對(duì)數(shù)運(yùn)算的原理如下所示：

lb x=k+lb （1+f）（3）

其中，對(duì)于lb （1+f）進(jìn)行分段線性逼近，絕對(duì)值模塊計(jì)算a_in的絕對(duì)值。前導(dǎo)零檢測(cè)模塊通過(guò)檢測(cè)前導(dǎo)零的個(gè)數(shù)來(lái)得到首1的位置k′。然后通過(guò)取首1之后的剩余位數(shù)來(lái)得到尾數(shù)值f。 f的高6位被用作尋址查找表來(lái)逼近lb （1+f）的非線性部分。其中，n表示尾數(shù)的位數(shù)，將首1的位置數(shù)減去尾數(shù)位數(shù)，將會(huì)得到首1的實(shí)際位置k。最終的轉(zhuǎn)換結(jié)果會(huì)通過(guò)連接k值和逼近區(qū)域部分的值來(lái)得到。

其中LUT存放查找表的系數(shù)，直接關(guān)系到最終消耗的硬件資源的多少。文獻(xiàn)[9]采用了15段的分段區(qū)間劃分，而經(jīng)過(guò)算法的自動(dòng)優(yōu)化，在同等精度4.1E-003條件下，本文的分段區(qū)間為7段，大大降低了查找表的資源消耗。

5結(jié)語(yǔ)

本文實(shí)現(xiàn)的改進(jìn)后的超越函數(shù)分段線性逼近算法可以在精度可控條件下，實(shí)現(xiàn)一個(gè)在不同取值范圍內(nèi)分段逼近區(qū)間的較優(yōu)劃分，能夠有效地降低分段區(qū)間的數(shù)量，減少查找表資源消耗。計(jì)算了多種不同的超越函數(shù)，并與其他文獻(xiàn)算法進(jìn)行對(duì)比，其結(jié)果充分表明：本文算法在保證計(jì)算精度的前提下，能夠有效地降低資源消耗，并且適用于不同的超越函數(shù)計(jì)算，而且其精度已知的特點(diǎn)也進(jìn)一步增加了算法的適用性和靈活性。同時(shí)，本文提出了這種改進(jìn)的算法，已經(jīng)用于自主開發(fā)的GPU上，并經(jīng)過(guò)FPGA驗(yàn)證。接下來(lái)的研究方向是采用查找表與更高階的多項(xiàng)式相結(jié)合的方式對(duì)超越函數(shù)進(jìn)行逼近，研究分段區(qū)間中的函數(shù)最優(yōu)逼近方式。

參考文獻(xiàn)：

[1]

焦繼業(yè)，穆榮，郝躍，等.面向移動(dòng)圖形頂點(diǎn)處理器的高性能低功耗定點(diǎn)特殊函數(shù)運(yùn)算單元[J].電子與信息學(xué)報(bào)，2011，33（11）：2764-2770.（JIAO J Y， MU R， HAO Y， et al. High performance and low power fixedpoint special function unit for mobile vertex processors [J]. Electronics & Information Technology， 2011， 33（11）： 2764-277.）

[2]

NICKLOLLS J， DALLY W J. The GPU computing era [J]. IEEE Micro， 2010， 30（2）： 56-69.

[3]

吳慶達(dá)，何書專，潘紅兵，等.32位浮點(diǎn)數(shù)正余弦函數(shù)FPGA實(shí)現(xiàn)方法[J].微電子學(xué)與計(jì)算機(jī)，2012，29（1）：113-116.（WU Q D， HE S Z， PAN H B， et al. Implementations of 32bits fixed and floating point trigonometric functions with FPGA [J]. Microelectronics & Computer， 2012， 29（1）： 113-116.）

[4]

鮑宜鵬.一種CORDIC算法優(yōu)化及32位浮點(diǎn)反正切函數(shù)FPGA實(shí)現(xiàn)[J].電子與封裝，2015，15（3）：22-25.（BAO Y P.One improved CORDIC algorithm of calculating 32bit floating the arctangent functions with FPGA [J]. Electronics and Packaging， 2015， 15（3）： 22-25.）

[5]

FLOREA C， VERTAN C. Piecewise linear approximation of logarithmic image processing models for dynamic range enhancement [J]. IEEE Transactions on Power Systems， 2011， 26（4）： 2581-2583.

[6]

王少軍，張啟榮，彭宇，等.超越函數(shù)FPGA計(jì)算的最佳等距分段線性逼近方法[J].儀器儀表學(xué)報(bào)，2014，35（6）：1209-1216.（WANG S J， ZHANG Q R， PENG Y， et al. Optimal equidistant piecewise linear approximation algorithm for the computation of transcendental functions in FPGA [J]. Chinese Journal of Scientific Instrument， 2014， 35（6）： 1209-1216.）

[7]

張建明，林亞平，傅明，等.傳感網(wǎng)絡(luò)中誤差有界的分段逼近數(shù)據(jù)壓縮算法[J].軟件學(xué)報(bào)，2011，22（9）：2149-2165.（ZHANG J M， LIN Y P， FU M， et al. Piecewise approximation based data compression algorithm with error bound in wireless sensor networks [J]. Journal of Software， 2011， 22（9）： 2149-2165.）

[8]

江釗.分段線性逼近法在梯級(jí)水電站優(yōu)化調(diào)度中的應(yīng)用[D].武漢：華中科技大學(xué)，2012：4-5.（JIANG Z. Application of piecewise linear approximation method to the optimal scheduling of casades hydroelectric plants [D]. Wuhan： Huazhong University of Science and Technology， 2012： 4-5.）

[9]

NAM B G， KIM H， YOO H J. Power and areaefficient unified computation of vector and elementary functions for handheld 3D graphics systems [J]. IEEE Transactions on Computers， 2008， 57（4）： 490-504.

[10]

DUNHAM J G. Optimum uniform piecewise linear approximation of planar curves [J]. IEEE Transactions on Pattern Analysis and Machine Intelligence， 1986， PAMI8（1）： 67-75.

[11]

BAJGER M， OMONDI A. Lowerror， highspeed approximation of the sigmoid function for large FPGA implementation [J]. Journal of Signal Processing Systems for Signal Image and Video Technology， 2008， 52（2）： 137-151.

[12]

周輝，李濤，邢啟江，等.數(shù)字曲線的線性逼近和分段識(shí)別[J].大連理工大學(xué)學(xué)報(bào)，1997，37（5）：576-580.（ZHOU H， LI T， XING Q J， et al. Linear approximation and piecewise identification of digital curves [J]. Journal of Dalian University of Technology， 1997， 37（5）： 576-580.）

[13]

NAM B G， KIM H， YOO H J. A lowpower unified arithmetic unit for programmable handheld 3D graphics systems [J]. IEEE Journal of SolidState Circuits， 2007， 42（8）： 1767-1778.