田征 杜慧敏 黄小康

摘要：针对超越函数计算中所采用的分段线性逼近算法存在的无法提前确定精度及部分区间资源浪费的问题，提出一种改进的分段线性逼近超越函数算法。该算法由预定义的逼近区间端点计算出用于逼近的线性函数，根据被逼近函数的凹凸性对所计算线性函数进行调整，在此基础上计算出预定义逼近区间内调整后函数与被逼近函数之间的最大误差；按照所需精度的要求，自动调整逼近区间，通过该过程的迭代，获得了较少分段次数。算法结果在Matlab上进行仿真，仿真结果表明，所提算法的分段数相比等分法减少了60%。所提算法在保证精度的前提下，降低了查找表（LUT）的资源消耗。

关键词：

分段线性逼近；超越函数；查找表；资源浪费；优化分段方法

中图分类号： TP391.75 文献标志码：A

0引言

图形处理器（Graphics Processing Unit， GPU）是各种嵌入式系统、个人机（Personal Computer， PC）、工作站和游戏机等不可缺少的重要部件。浮点超越函数单元是GPU数据通路中的重要部件，其性能直接影响图形渲染效果[1-2]。

目前在现场可编程门阵列（FieldProgrammable Gate Array， FPGA）中计算超越函数常用的方法有级数近似法、查表法（Look Up Table， LUT）、坐标旋转数字计算（Coordinate Rotation Digital Computer， CORDIC）算法和分段线性逼近法等。

其中，级数近似法展开式较为复杂，硬件实现复杂，资源消耗较大。查表法虽然计算简单、易于实现，但是所需存储单元随着计算精度的提高呈指数形式增加，资源消耗大[3]。CORDIC算法[4]作为一种便于FPGA实现的超越函数计算方法得到了广泛关注，但其收敛速度与数据的表示精度成反比，当精度要求较高时，算法的迭代次数较多，计算延迟会增大。与之相比，分段线性逼近法[5-6]将低阶多项式与较小的查找表相结合，资源消耗少，速度较快，被广泛应用于传感网络中的数据压缩[7]、非线性模型到线性模型的转换[8]以及图形图像处理[9]等领域，已成为在计算资源有限条件下超越函数计算的一种较为理想的选择。

为了提高函数的计算速度，在硬件设计中一般利用查找表实现分段线性逼近算法[6]，段数越多，查找表越大，误差则越小，因此，分段线性逼近算法需要在查找表、精度和分段数目之间寻求一个合理的平衡。其中有均方误差法[10]、区间2k等分法[11]、面积法[12]等多种分段算法，但这些算法都是在计算完成后才可以知道分段精度，难以在精度已知的条件下完成分段数的计算。

文献[6]中提出了一种超越函数计算的最佳等距分段线性逼近（Optimal EquiDistant PieceWise Linear approximation，OED_PWL）方法，该方法可以通过调节分段数目来完成计算精度的灵活控制，从而可以在保证精度的条件下完成分段数的计算。相较区间2k等分法其算法性能、资源消耗等方面的表现更好；但是，根据文献[6]提出的方法进行tanh函数计算时，部分区间的计算误差较大，原因是tanh函数斜率在某些分段较平缓，而在另外一些分段斜率较陡峭，若采用均匀分段的方法，则存在分段平缓处的资源浪费问题。在后面的实验结果中，本文会与OED_PWL算法进行对比。

因此本文提出一种改进的超越函数分段线性逼近方法，通过自动识别逼近误差进行分段使其不断细分，实现了一个精度已知条件下较优的分段方法。

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

分段线性逼近算法的原理是把非线性特性曲线分成若干个区段，在每个区段中用直线段近似地逼近特性曲线。算法原理如图1所示。

分段线性逼近算法是根据一定的分段方法，如面积法、等分法等，将曲线f（x）所在的一段区间[A，B]划分为若干段，在每一段里面通过一个线性函数来进行逼近：

y（x）ax+b（1）

其中：a表示线性函数的斜率；b表示线性函数的偏移量。其中，在每一段逼近区间[x1，x2]里面，区间的两个端点需要满足：

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

在通过某种方法得到[A，B]上超越函数曲线的分段情况以及各个段的逼近系数ai、bi后，即可完成超越函数的分段线性逼近。

算法硬件实现分为3步：

步骤1按照某种方法进行分段区间的划分，求取所有分段点。

步骤2在每个分段区间，根据两个端点的值求取直线段f（x）=ax+b，将每一段的ai和bi存放在查找表中以待使用。

步骤3利用查找表，读取每一段直线方程的系数，通过乘法和加法计算每一段函数的逼近值，从而实现超越函数的分段线性逼近。

2改进的分段线性逼近算法

2.1算法原理

如何在给定的精度下，得到一个较好的分段，使得应用分段计算值满足给定的精度要求是分段线性逼近算法的一个难题。本文提出一种基于区间2k等分法，通过迭代来获得分段区间的分段逼近方法，接下来对其原理进行详细介绍。

区间2k等分法会将[A，B]上的曲线均匀地分为2k段，k的取值取决于查找表的大小及资源的分配情况。这是分段线性逼近中最常用的一种分段方法，具有实现简单、精度较高等优点；但是由于超越函数曲线的斜率在定义域内是不均匀分布的，部分计算域内曲线段斜率较大，因此采用直线段线性逼近误差较大；而另外一部分计算域曲线段斜率较小，逼近误差较小。为了满足误差较大曲线段的误差要求，往往会导致分段数过多，导致斜率较小曲线段的资源浪费。

因此本文提出一种改进型的超越函数分段线性逼近算法，可以在精度已知的前提下，完成对分段数的最优选取，使得计算函数误差在一定的范围内。算法的主要思想为：对斜率较大、误差较大的函数区间进行进一步细分，避免了对曲线平滑部分分段的浪费，减少了查找表资源的消耗。

以对数lb x在区间[1，5]的计算说明算法的思想：首先整个逼近区间等分为[1，3]和[3，5]回复：分段区间没有问题两段，并对两段分别进行线性逼近，得到line1的两条直线段，可以看出，在[1，3]区间的直线逼近误差比较大，因此，对区间[1，3]继续进行细分，分为更小的两个区间：[1，2]和[2，3]在等于2时，两个区间都包含在内，其中一个需用开区间来表示，请明确。回复：分段区间没有问题，并对其分别进行线性逼近，得到两条误差较低的直线段line2，完成分段区间的细分过程。

2.2算法实现流程

本文实现的算法流程如图3所示，通过迭代的方法来完成分段区间的计算。

图4说明了算法的一次逼近过程，其中逼近区间段的中点处逼近误差为d，假设原逼近直线y0=ax+b向曲线方向偏移d/2，得到新的逼近直线y1=ax+b+d/2，新的直线更加逼近曲线。

本文在求取被逼近函数f（x）与新的逼近直线y1 （x）之间的最大误差dmax=f（x）-y1（x），对dmax求导数，令导数等于零，将得到x的值，代入dmax计算出最大误差。

在本文的改进型分段线性逼近算法中，所有数据采用双精度类型，在每个分段区间上，判断被逼近函数的凹凸性以及该函数在区间内是否存在拐点，求取逼近直线以及被逼近函数之间的最大误差，然后与给定误差进行比较：若小于给定误差，则停止该段的继续细分；若大于给定误差，则对该段逼近区间继续细分。

若删除，可从此处删除伪代码。

下面是改进型的分段线性逼近算法的C语言伪代码。

图3的表述，与此处的“c语言伪代码”？是否有重复表达现象，是否可以删除某一种表述方式？请明确。

回复：其中分段区间没有问题，还有就是那个关于伪代码和流程图重复问题，我和我们老师商量了，可以删除也可以不删除，如果要删除，就把伪代码部分删除。（附件是删除伪代码的稿件）

：

程序前

void subsection（double x0， double x1）

{

double y0=x0坐标处要逼近的函数值f（x0）；

double y1=x1坐标处要逼近的函数值f（x1）；

y′=f ″（x）=sign/*对逼近函数求导判断凹凸性，sign=1表示凸函数，sign=0表示凹函数*/

令f ′（c）=0求出c值，c则为函数逼近区间内的拐点；

if（y′=0）/*函数有拐点*/

{

重新划分逼近区间为[x0，c]及[c，x1]，确保每段新的逼近区间都是单调的；

}

else；

mid=（x0+x1）/2.0；

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

b=y1-a*（x1）；

max=d=f（mid）-（a*mid+b）；/*mid坐标处的逼近误差*/

if（sign==1）

b1=b+d/2.0；

else

b1=b-d/2.0；

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

d=f（x）-y//被逼近函数与逼近函数的差值

d′=f ′（x）-y′//对函数d求导数

d′=f（m）-y1′（m）=0//求函数的极值点m

max=f（m）-y1（m）；//逼近区间内的最大误差

if （max<指定误差）{

打印x0，x1坐标，并打印逼近直线参数m及n到指定文档；

}

else {

subsection（x0，mid）；

subsection（mid，x1）；

}

}

//说明：x0与x1为要逼近的函数取值范围；

程序后

3实验结果与对比

为了更好地说明问题，采用1E-3和1E-4作为预设的精度，可以满足移动端图形图像处理等应用的精度要求；取值范围设定为[1，2]，可以充分地说明逼近的正确性。表1为各个超越函数在区间2k等分法和本文算法下得到的分段数量对比。

可以看出，改进后的分段线性逼近算法对斜率变化较大的超越函数作用会比较明显。而且如果精度提升，2k分段数会剧烈增长，而本文提出的算法分段数增加较慢。

表2为本文算法与其他两种算法在相同实验条件下的分段结果比较，可以看出：在保证误差精度的前提下，本文算法的分段数相比等分法减少了60%以上，相比OED_PWL法减少了20%。对于FPGA硬件实现而言，分段数越少，需要的查找表则越小，从而节省了硬件查找表资源。

针对超越函数中的对数函数lb （1+x），文献[9，13-14]等多篇文献都对其进行了分段区间的优化。表3为本文算法对对数函数的分段逼近结果与其他文献算法的对比，由于其他文献的分段算法没有做到精度已知下的分段划分，因此实现精度各不相同。

由表2与表3的实验结果可以得出：

1）与其他已经提出的算法相比，在相同精度下，本文的算法划分的逼近区间更少，则查找到所属分段的时间就会变少，所占据的存储空间更小，更适合用硬件实现。

2）根据表3，与文献[9]相比，分段数目相差不多的情况下，本文的算法可以大大降低逼近误差。

综上所述，改进型的分段线性逼近算法同样适用于tanh（x）函数、lb x函数等多种超越函数算法。

4对数函数的FPGA实现

完成算法的基本设计之后，本文将其应用在某个超越函数运算单元中，用于计算以2为底的对数转换结果，其电路结构如图5所示。

对数运算的原理如下所示：

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

其中，对于lb （1+f）进行分段线性逼近，绝对值模块计算a_in的绝对值。前导零检测模块通过检测前导零的个数来得到首1的位置k′。然后通过取首1之后的剩余位数来得到尾数值f。 f的高6位被用作寻址查找表来逼近lb （1+f）的非线性部分。其中，n表示尾数的位数，将首1的位置数减去尾数位数，将会得到首1的实际位置k。最终的转换结果会通过连接k值和逼近区域部分的值来得到。

其中LUT存放查找表的系数，直接关系到最终消耗的硬件资源的多少。文献[9]采用了15段的分段区间划分，而经过算法的自动优化，在同等精度4.1E-003条件下，本文的分段区间为7段，大大降低了查找表的资源消耗。

5结语

本文实现的改进后的超越函数分段线性逼近算法可以在精度可控条件下，实现一个在不同取值范围内分段逼近区间的较优划分，能够有效地降低分段区间的数量，减少查找表资源消耗。计算了多种不同的超越函数，并与其他文献算法进行对比，其结果充分表明：本文算法在保证计算精度的前提下，能够有效地降低资源消耗，并且适用于不同的超越函数计算，而且其精度已知的特点也进一步增加了算法的适用性和灵活性。同时，本文提出了这种改进的算法，已经用于自主开发的GPU上，并经过FPGA验证。接下来的研究方向是采用查找表与更高阶的多项式相结合的方式对超越函数进行逼近，研究分段区间中的函数最优逼近方式。

参考文献：

[1]

焦继业，穆荣，郝跃，等.面向移动图形顶点处理器的高性能低功耗定点特殊函数运算单元[J].电子与信息学报，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]

吴庆达，何书专，潘红兵，等.32位浮点数正余弦函数FPGA实现方法[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算法优化及32位浮点反正切函数FPGA实现[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]

王少军，张启荣，彭宇，等.超越函数FPGA计算的最佳等距分段线性逼近方法[J].仪器仪表学报，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]

张建明，林亚平，傅明，等.传感网络中误差有界的分段逼近数据压缩算法[J].软件学报，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]

江钊.分段线性逼近法在梯级水电站优化调度中的应用[D].武汉：华中科技大学，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]

周辉，李涛，邢启江，等.数字曲线的线性逼近和分段识别[J].大连理工大学学报，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.