位宽优化中乘法运算的一种自动范围分析方法

孙瑞一，张 岩

(哈尔滨工业大学深圳研究生院网络环境智能计算重点实验室，518055 广东深圳)

位宽优化中乘法运算的一种自动范围分析方法

孙瑞一，张 岩

(哈尔滨工业大学深圳研究生院网络环境智能计算重点实验室，518055 广东深圳)

乘法是硬件平台中最基本的非线性运算，而且在自动位宽优化过程中，目前的范围分析方法没有在精确的范围分析结果和计算复杂度之间做很好折衷.为了在较低的计算复杂度前提下更准确地分析乘法运算结果的范围，提出了改进的仿射近似法(NAA).在改进的仿射近似法中，利用额外噪声项来表示近似产生的误差，并根据误差的特点把误差分成两部分，在不增加计算复杂度的前提下更准确地估计误差的范围.新方法的计算复杂度是O(M1)，其中M1是乘法的两个操作数中非零噪声个数的和.实例分析表明，利用该方法得到的乘法结果范围的准确程度是用简单估计法得到的准确程度的1.47倍，和切比雪夫近似法的准确度接近.

位宽优化;范围分析;乘法;仿射算术;仿射近似法

在硬件系统的设计中，为了得到较高的数据精度及动态范围，算法模型都用浮点数表示.但是在硬件上实现浮点运算的代价很大，为了提高运行速度、降低功耗、节省面积，在硬件实现阶段数据一般都采用定点数表示.将浮点数转化为定点数的过程被称为位宽优化，也叫做定点化.位宽优化的目的是在满足系统规格要求的前提下寻找最佳的定点数位宽组合使系统的代价(面积、功耗、速度等)最小.位宽优化包括范围分析(对整数部分的优化)和精度分析(对小数部分的优化)两个步骤.

位宽优化是NP－hard问题［1］，常用的分类方法有动态法［2－6］和静态法［7－10］.动态法通过对大量的测试向量进行反复的蒙特卡洛仿真来确定信号的位宽.它的结果质量较高，但优化时间长，甚至占到整个设计周期的50%以上［11］.而且，动态法不能保证测试向量没有覆盖到的情况也满足系统规格要求.静态法使用代码分析手段来推导出信号的位宽，它不需要测试向量，所以它的执行速度快、人为干扰因素小，而且简便易行，适于大规模系统的自动优化分析.

在硬件平台，如Xilinx公司的FPGA或全定制的ASIC芯片中，乘法运算是最基本也是最常用的非线性运算.其他的非线性运算，如除法、余弦、对数等超越函数一般都是通过加、减、乘这些基本运算计算的.但是目前乘法运算的范围分析的静态法得不到范围的精确解，所以对乘法运算结果的位宽优化的准确程度就关系到整个电路的位宽优化的准确程度.本文讨论基于静态法的对乘法结果的范围分析方法.

区间算术(interval arithmetic，IA)是1960年由Moore［12］提出的解决范围分析的方法.Cmar等［13］首次采用IA对信号的变化范围进行分析.但是IA在估计信号范围时过于保守，甚至不切实际.

仿射算术(affine arithmetic，AA)保留了区间之间的相关性，这使得它比IA更精确.AA非常适于对线性操作结果的范围分析.但是AA不能为乘法等非线性操作提供准确的仿射形式.为解决该问题，Stolfi等［14］提出了乘法的仿射近似法，包括简单估计法和切比雪夫近似法.简单估计法计算效率高，但误差较大，分析范围最多可达到真实范围的4倍.这种误差沿着数据流累积到输出，会导致误差爆炸，这限制了它在大系统中的应用;切比雪夫近似法利用比简单估计法更精确的仿射形式，实现较好的范围分析结果，但是它的计算过程太复杂所以不适合在大系统中应用.

为了解决计算效率和计算复杂度之间的问题，Zhang等［15］提出了一种叫N-级近似的，在分析精度和复杂度之间作折中的方法.其分析结果比用简单估计法得到的准确，比切比雪夫近似法简单.但是仍不能满足大系统对分析范围的精度和复杂度的要求.Pang等［16］提出了一种新的范围分析方法，它混合了传统的IA、AA和算术变换(arithmetic transform，AT)，能得到比 AA更精确的范围结果，但是执行时间比AA长很多.所有这些对AA的改进方法，都以牺牲计算复杂度为代价来换取分析的精度.它们都是把自变量的仿射形式看成是一个整体，没有考虑仿射形式中的每一个噪声，这样就忽略了不同噪声之间的独立性和相同噪声在不同变量中的相关性.

为了更简便、更精确的实现对乘法运算结果的范围分析，本文提出一种名为改进的仿射近似法(novel affine arithmetic approximation，NAA)的新的仿射近似法.NAA利用一个仿射形式近似的表示乘法运算的结果，并利用增加的额外噪声项来表示近似产生的误差.为了利用现有方法忽略的独立性和相关性，这个误差用各个噪声来表示，并且在不增加计算复杂度的前提下更准确的估计误差的范围，从而减少了新增噪声的系数.这种方法能比简单估计法得到更准确的分析范围，而且计算复杂度和简单估计法一样，比切比雪夫近似法更简单.

1 仿射算术

为了更清楚的分析仿射近似的过程和更方便的推导NAA，首先介绍仿射算术的基础知识.

仿射算术用一次多项式的仿射形式^x表示一个未知的信号x为

式中:εi=［-1，1］.对于未知信号 x，x0表示它的中心值;εi为第i个噪声，代表信号x的一个独立的、未知的噪声源;xi为这个噪声的系数;xiεi、x0+xiεi分别为噪声项.

令x0=(xmax+xmin)/2，x1=(xmax-xmin)/2，信号的范围区间¯x=［xmin，xmax］可以转化为等价的仿射形式为

AA通过与计算链中的其他信号共享噪声εi来保留信号之间的相关性.

乘法的简单估计法的仿射形式为

这种方法的分析结果误差很大，但是其计算复杂度是O(M1)，其中M1=max(n1，n2)，n1、n2分别为信号x、y中非零噪声的个数.

式中:m、n分别为(x-x0)(y-y0)的最大值和最小值.

计算最值m和n的计算复杂度是O(M2log M2)，其中 M2=n1+n2.

2 乘法的改进的仿射近似法

2.1 算法分析

自变量为仿射形式的乘法，可以表示为

因为z不是仿射形式，选择一个仿射形式fz来近似的表示它.式(1)的前3项组成了一个仿射形式，为了简单，把这个仿射形式作为近似的仿射形式，即

式中:fz为z的近似表示，所以fz与z之间存在误差，引入一个新的独立噪声项来表示这个误差df，即

计算df的真实最大值和最小值的计算量难以承受，所以dfmax、dfmin分别为df的近似最大值和近似最小值.dfmax、dfmin与真实的最大值和最小值越接近，^z的形式越准确.

综上所述，z的仿射形式可以记为

2.2 仿射形式

df可以表示为

式中:εi，εj=［-1，1］.估计df的近似的最大值和最小值的最简单的方法就是把εi，εj的最大值和最小值代入.根据式(2)，当 i=j时，有 xiyjεiεj=.把εi，εj=［-1，1］代入，xiyi估计的近似的最大值和最小值比xiyjεiεj更精确，而且没有增加计算复杂度.因此把df分成两部分，即

距离df可以记为

d1的最大值和最小值分别为

和它等价的仿射形式为

d2的最大值和最小值分别为

和它等价的仿射形式为

计算一次乘法要引入两个独立噪声会使接下来的计算噪声越来越多.所以需要化简根据式(4)，的范围为

和它等价的仿射形式为

式(4)中的εn+1和εn+2之间有相关性，但是将df分为d1和d2并没有考虑它们之间的相关性，所以df的真实范围一定比的范围小，这是过估计.因此，NAA的准确程度比切比雪夫近似法的准确程度稍低.

2.3 乘法的改进的仿射近似法的表达式及计算复杂度分析

有了式(2)和新引入的式(5)，乘法的改进的仿射近似法的表达式为

设n1、n2分别为^x、^y中非零噪声的个数，M1=max(n1，n2)为n1和n2的最大值，M2=n1+n2，则切比雪夫近似法的计算复杂度为 O(M2log M2)［14］.简单估计法的计算复杂度为O(M1).

为了计算NAA的计算复杂度，把式(6)改写为

式(7)的计算复杂度为

可以看出，NAA的计算复杂度和简单估计法的一样，都是O(M1).很明显M1＜M2，所以有O(M1)＜O(M2log M2).根据计算复杂度的分析，NAA比切比雪夫近似法要简单.

3 实例分析

为了比较NAA的性能，本文以包含乘法的计算为例.第1个实例是多项式近似.第2个实例是B-splines.这两个实例均来自于文献［8］.第3个实例是多变量多项式，来自于文献［17－19］.将仿射算术用C++实现，乘法分别采用简单估计法、切比雪夫近似法和改进的仿射近似法.这些实例的计算平台是内存为2.99 G的Intel(R)Pentium(R)CPU G630@2.7 GHz的32位双核PC.

3.1 实例

实例1 多项式近似

对 y=ln(1+x)，x=［0，1］的多项式近似进行信号范围分析，采用文献［20］的4级多项式表达为

其中，5个四舍五入到小数点后第4位的系数由多项式曲线拟合的方法获得.

实例2 B-splines

B-splines常用于图像卷绕，它的基本函数B0、B1、B2、B3分别定义如下，其中，输入 u= ［0，1］.

实例3 变量多项式

考察表1中的8个多变量多项式.

表1 多变量多项式函数及输入范围

3.2 3种方法分析范围的比较

利用仿射近似法计算得到的变量范围用分析范围来表示.用分析范围比率来表示分析范围的准确程度.分析范围比率是用3种方法得到的分析范围的区间间隔(ymax-ymin)与真实范围的区间间隔的比值，表示了分析范围的区间间隔是真实范围的区间间隔的倍数.其中，ymax、ymin分别为变量的最大值和最小值.分析范围比率越接近1说明分析范围越接近真实范围.

表2给出了用3种方法计算得到的3个例子的分析范围和范围比率.它显示了利用NAA计算得到的分析范围全部覆盖了真实的输出范围.对于这些例子，利用简单估计法计算得到的范围比率从1.04～281.19;利用切比雪夫近似法计算得到的范围比率从1.03～233.66;利用NAA计算得到的范围比率从1.03～221.78.利用NAA得到的范围比率是利用简单估计法得到的范围比率的0.18～0.99倍，范围比率倍数的平均是0.68倍;利用NAA计算得到的范围的比率是利用切比雪夫近似法计算得到的范围比率的0.33～1.39倍，范围比率倍数的平均是0.99倍.这说明，平均来看，利用NAA计算得到的分析范围的准确程度是用简单估计法计算得到的分析范围的准确程度的1.47倍，和切比雪夫近似法的分析范围的准确度接近.

表2 3种方法分析范围和范围比率的比较

对于实例3中的第7个函数，利用3种静态法得到的范围比率都在200以上.这是因为这个函数最终的乘法中的两个乘数都包含很多乘法，每一次乘法运算都会产生误差，而且这两个乘数都是相同自变量的函数，它们的相关性很强，这种情况下两个乘数的误差会相互放大，使最终得到的范围结果比真实范围大很多.所以静态法不适用于对这种情况的范围分析，动态法会得到更准确的分析范围.

从表2可以看出，这3种方法有一个共同的不足之处，即分析范围的下限误差较大.这是因为当变量x用仿射形式表示时，x的范围是相对于x0中心对称的.在计算过程中，x0是由计算链中变量决定的，很难保证它是在x范围的中心.当^x要满足变量范围的上限或下限时，由于x0不在x范围的中心，^x表示的x范围的下限或上限就有了较大的误差.在本文所举的例子中，都是上限大于下限的，所以表现为分析范围下限误差较大.

3.3 3种方法整数位宽和执行时间的比较

变量的整数位宽可以由变量范围推导.表3列出了3种方法得到的输出变量位宽和CPU运行时间.B-splines的4个基本函数是在同一个程序中完成的，所以只有一个CPU时间.NAA最多比简单估计法节省2个bit位宽，比切比雪夫近似法最多节省1个bit位宽.与简单估计法相比，NAA需要0.95～1.39倍的CPU时间，而切比雪夫近似法需要54.01～57.97倍的 CPU时间.NAA和简单估计法有着近似的执行时间，不到切比雪夫近似法执行时间的1/50.

表3 3种方法得到的整数位宽和CPU时间的比较

4 结论

1)利用一个仿射形式近似的表示乘法运算的结果，并利用额外噪声项来表示近似产生的误差.

2)把近似产生的误差分成两部分，分别用两个噪声来表示，能得到更精确的近似误差的范围.

3)与以往的范围分析的改进方法相比，新方法在不增加计算复杂度的前提下的，更准确的估计变量的范围.实例分析显示NAA的分析范围的准确程度和切比雪夫近似法的准确程度接近，比简单估计法更准确.而且，它的计算复杂度和简单估计法的一样，都是O(M1)，比切比雪夫近似法的小.

［1］CONSTANTINIDES G，WOEGINGER G.The complxity ofmultiple wordlength assignment［J］. Applied Mathematics Letters，2002，15(2):137－140.

［2］KUM K I，SUNG W.Combined word-length optimization and highlevel synthesis of digital signal processing systems［J］.IEEE Trans on Computer-Aided Design Integrated Circuits，and Systems，2001，20(8):921－930.

［3］CAFFARENA G，CARRERAS C，LOPEZ J A，et al.SQNR estimation of fixed-point DSP algorithms［J］.Eurasip Journal on Advances in Signal Processing，2010，2010(21):1－12.

［4］朱珂，华林，周晓方，等.JPEG2000中小波滤波器的定点分析及其VLSI实现［J］.固体电子学研究与进展，2004，24(4):466－471，504.

［5］马志强，季振洲，胡铭曾.基于超窄数据的低功耗数据Cache方案［J］.计算机研究与发展，2007，44(5):775－781.

［6］BANCIU A，CASSEAU E，MENARD D，et al.Stochastic modelingforfloating-pointto fixed-point conversion［C］//Proceedings of IEEE Workshop on Signal Processing Systems(SiPS).Beirut，Lebanon:IEEE Computer Society Press，2011:180－185.

［7］SARBISHEI O，RADECKA K，ZILIC Z.Analytical optimization of bit-widths in fixed-point LTI systems［J］.IEEE Trans on Computer-Aided Design of Integrated Circuits and Systems，2012，31(3):343－355.

［8］LEE D－U，GAFFAR A A，CHEUNG R C，et al.Accuracy guaranteed bit-width optimisation［J］.IEEE Trans on Computer-Aided Design of Integrated Circuits and Systems，2006，25(10):1990－2000.

［9］ROCHER R，MENARD D，SCALART P.Analytical approach for numerical accuracy estimation of fixed-point systems based on smooth operations［J］.IEEE Trans on Circuits and Systems I-Regular Papers，2012，59(10):2326－2339.

［10］KINSMAN A B，NICOLICI N.Computational vectormagnitude-based range determination for scientific abstract data types［J］.IEEE Trans on Computers，2011，60(11):1652－1663.

［11］KEDING H，WILLEMS M， COORS M， et al.FRIDGE:a fixed-point design and simulation environment［C］//Proceedings of Design，Automation and Test in Europe.Paris:DATE，1998:429－435.

［12］MOORE R E.Interval Arithmetic and Automatic Error Analysis in Digital Computing［D］.California:Stanford University，1962.

［13］CMAR R，RIJNDERS L，SCHAUMONT P，et al.A methodology and design environment for DSP ASIC fixed point refinement［C］//Proceedings of Design，Automation and Test in Europe.Munich:ACM，1999:271－276.

［14］STOLFI J，de FIGUEIREDO L H.Self-validated Numerical Methods and Applications［M］.Rio De Janeiro:Brazilian Mathematics Colloquium monograph，IMPA，1997.

［15］ZHANG Linsheng，ZHANG Yan，ZHOU Wenbiao.Tradeoff between Approximation accuracy and complexity for range analysis using affine arithmetic［J］.Journal of Signal Processing Systems，2010，61(3):279－291.

［16］PANG Y，RADECKA K.An efficient algorithm of performing range analysisforfixed-pointarithmetic circuits based on SAT checking［C］//Proceedings of IEEE International Symposium on Circuits and Systems(ISCAS).Rio de Janeiro，BRAZIL:IEEE Computer Society Press，2011:1736－1739.

［17］SHEKHAR N，KALLA P，ENESCU F.Equivalence verification of arithmetic datapaths with multiple wordlength operands ［C］//Proceedings of Design，Automation and Test in Europe.Munich:DATE，2006:824－829.

［18］GOPALAKRISHNAN S，KALLA P，MEREDITH M B，et al.Finding linear building-blocks for RTL synthesis of polynomial datapaths with fixed-size bit-vectors［C］//Proceedings of International Conference on Computer-Aided Design.San Jose:IEEE Computer Society Press，2007:143－148.

［19］SHOU Huahao，SONG Wenhao，SHEN Jie，et al.A recursive taylor method for ray casting algebraic surfaces［C］//Proceedings ofInternationalConference on Computer Graphics and Virtual Reality，Las Vegas:CGVR，2006:196－204.

［20］HORNER W G.A new method of solving numerical equations of all orders，by continuous approximation［J］.Philosophical Transactions of the Royal Society of London，1819，109:308－335.

A range analysis in automatic word length optimization for multiplication

SUN Ruiyi，ZHANG Yan
(Key Laboratory of Network Oriented Intelligent Computation，Shenzhen Graduate School，Harbin Institute of Technology，518055 Shenzhen，Guangdong，China)

To achieve more accurate result and lower computational complexity of range analysis for multiplication in automatic word length optimization，this paper presents a novel refined affine approximation method of multiplication for range analysis in automatic word length optimization，which is named novel affine arithmetic approximation(NAA).In NAA，a new noise term represents the error which is caused by approximation.This error is estimated more accurately without increasing the computational complexity.The computational complexity of NAA is O(M1)，where M1denotes the total of the nonzero noise of the two multipliers.In experiments，the accuracy of the range using NAA is 1.47 times of that using trivial range estimation，and the same as that using Chebyshev approximation.

word-length optimization;range analysis;multiplication;affine arithmetic;affine approximation method

TP391.72

A

0367－6234(2014)03－0043－06

2013－02－27.

深圳市科技研发基础研究计划资助项目(JC201005260168A).

孙瑞一(1980—)，女，博士研究生;

张 岩(1969—)，男，教授，博士生导师.

张 岩，ianzh@foxmail.com.

(编辑 张 红)