基于CORDIC的反双曲正切函数的FPGA实现

1、 CN4321258/TP ISSN 10072130X 计算机工程与科学COMPU TER EN GIN EERIN G &SCIENCE2009年第 31卷第 5期 Vol 131,No 15,2009文章编号 :10072130X (2009 0520150203基于 CO RD IC 的反双曲正切函数的 FP GA 实现 3 The FP GA Implementatio n ofHyperbolic Arctangent Functio n陈石平 , 李 全 , 莫丽兰 , 段吉海CHEN Shi 2ping ,L I Q u an ,MO Li 2lan ,DUAN Ji

2、2hai(桂林电子科技大学信息与通信学院 , 广西 桂林 541004(School of Inform ation and Communications , G uilin U niversity of E lectronic T echnology , G uilin 541004,China摘 要 :双曲函数的应用领域十分广泛 。本文首先介绍 CORDIC 算法双曲系统的基本原理及其计算模式 , 对 CORDIC 内核及其处理单元做了详细分析 。 在迭代算法的基础之上 , 采用流水线技术 , 以面积换速度 , 给出了一种基于流 水线的 CORDIC 来实现反双曲正切函数 , 具有很高的精

3、度和很快的速度 , 使设计出的软核能够在精度要求很高的场合中 运行 。 用 Verilog HDL 对其编程设计和进行功能仿真 、 时序仿真及下载测试的结果表明 , 该函数具有很好的实用性 。 Abstract :Hyperbolic f unctions are widely applied in many areas. This paper first the principle of the CORDIC algorithm hyperbolic system and its computation model , and gives of pipelining the CORDIC k

4、ernel and its processor unit for the design. Based on the , implementing hyper 2 bolic f unctions based on pipelining CORDIC is proposed. It speed at the cost of FP GA resources , making the soft core adapted to 22Finally the results of programming the CORDIC , performing the f ,and downloading it t

5、o the FP GA show that it features great practicability.关键词 :; ; 流水线K ey w ords :CORDIC algorithm ;hyperbolic arctangent ;pipelining中图分类号 :TP322. 2文献标识码 :A1 引言双曲函数的应用领域十分广泛 , 涉及图像信号处理 、 滤 波技术 、 FFT 变换 、 LDPC 码软判决译码等领域 。因此 , 设 计并实现双曲函数显得尤为重要 , 用硬件实现超越函数的 算法 , 按照数学公式和对应的实现方式的不同 , 可以分为查 表法 、 多项式近似法 、 基于

6、查表的多项式结合方法 、 有理数 近似和逐位法五类 。经过对这些算法进行分析和比较权 衡 , 本文采用 CORDIC 算法来实现双曲函数中的反双曲正 切函数 , 并用 FP GA 实现硬件测试 。CORDIC (Coordinate Rotation Digital Computer , 简称 CORDIC 算法即坐标旋转数字计算方法 , 是 Volder J D 1于 1959年首次提出的 , 主要用于三角函数 、 双曲函数 、 指数 函数 、 对数函数的计算 。 该算法通过基本的加和移位运算 代替乘法运算 , 使得矢量的旋转和定向计算不再需要三角 函数 、 乘法 、 开方 、 反三角 、

7、反双曲正切等函数 。2 COR DIC 算法原理2. 1 双曲系统CORDIC 算法包含圆周系统 、 线性系统 、 双曲系统三种 旋转系统 2。 本文采用双曲系统 , 其推导如下 :由于 e =sinh +cosh , e 只需计算 sinh 和 cosh 即 可 , 又由于 :cosh (± =cosh cosh ±sinh sinh sinh (± =sinh cosh ±cosh sinh (1 图 1为平面坐标旋转 , X i =cosh , Y i =sinh , 从中051 3收稿日期 :2008202201; 修订日期 :200820520

8、7基金项目 :广西研究生教育创新项目 (2007105950810m15作者简介 :陈石平 (19812 , 男 , 湖南郴州人 , 硕士生 , 研究方向为 FP GA 和无线通信 ; 段吉海 , 硕士 , 副教授 , 研究方向为通信技术、 EDA 和专用集成电路设计。Address :Unit 322,10Fuhouli Lane , Guilin , Guangxi 541001,P. R. China 图 1 坐标旋转图可以看出 , 将向量 (X i , Y i 的角度 旋转 角度 , 得到一 个新的向量 (X j , Y j , 那么有 :X j =cosh (+ =cosh cosh

9、 +sinh sinh= X i cosh +Y i sinh Y j =sinh (+ =sinh cosh +cosh sinh = X i sinh+Y i cosh(2 将式 (2 化为矩阵表达形式 , 平面旋转定义如下 :X j Y j=coshsinh sinh coshX i Y i (3 使用迭代的方法 , 旋转角度可在多步完成 , 每一步的旋 转只完成其中的一小部分 , N 步旋转将会完成一个 角度旋转 。 由式 (3 得到单第 n 步旋转以及提取 cosh n 因子 , 有式 (4: X n+1Y n+1=cosh n sinhnsinh ncosh X n Y = cos

10、h n1tanhntanh n1 X n Y (4 由 于 cosh=2/(cosh 2-sinh 2=1/-tanh 2, 式 (4 中 每 一 步 的 旋 转 角 度 ntanh -12-n。 由于 tanh -1(20 , 开始 , 而且 n =1S n n=,±, 角度 。 S n , 有 tanh n =S n 2-n, 即 :X n+1Y n+1=cosh n1S n 2-S n 2-n1X n Y (5 再通过 N 次迭代后 :Xj Y j=Nn =1cosh n1S n 2-S n 2-n1X i Y i= K 3Nn =11S n 2-S n 2-n1X i Y (

11、6 将比例因子从迭代公式中提取出来 , 定义 K 为增益因 子 , 则有 :K =1/P =Nn =1cosh tanh -12-n=n =1-2-2n 0. 8297816(7 max =j =1tanh -1(2-j 0. 9578885(8 具体增益 K 取决于迭代次数 N 。对于所有的初始向量和旋转角度而言 , K 是一个常数 , 通常把 K 称作聚焦常 数 2。 式 (7 中常数 P 为 K 的倒数 , 式 (8 中 max 为双曲系 统中不重复迭代的最大旋转角度 。2. 2 计算模式双曲系统有旋转模式和向量模式 3两种计算模式 。在旋转模式中 , Z 为初始化需要旋转的角度 , 当

12、 Z 旋转变为 0时 , 用来计算双曲正余弦函数 。 本文采用向量模式 :它跟旋 转模式不同 , 其区别是旋转的方向取决于 Y 而不是 Z 的符 号 , 即将输入向量 Y 变为 0, 计算相应的 tanh -1(反双曲正 切 值 。 Y 为初始化值 , 当 Y 旋转变为 0时 , 式 (5 的输出变 为 :X n+1=X n +S n Y n 2-n Y n+1=Y n +S n X n 2-nZ n+1=Z n -S n tanh -1(2-n (9式 (9 中 , Y n <0时 , S n =+1, 反之 S n =-1。 经过 N 次 迭代后 ,CORDIC 的输出变为 3:X

13、n+1=K 02-Y 02Y n+1=0Z n+1=Z 0+tanh -1(X 0(10其中 , P =Nn =1-2-2n。如果 Z 0=0, 对 于 给 定 的 X 0和 Y 0, N 次 迭代 后CORDIC 公式的输出变为 :X n+1, Y n+1, Z n+=K02-Y 02, 0, tanh -1(X 0(11 分析式 (11 可得 。 对式 (10 进行函 x 0y , 可以得到不同 , , 如表 1所 , R 和 z R 为函数输出值 。表 1 其它的函数计算 模式 x 0 y 0 z 0 x R z R 备注向量1a0-a 2z R =tanh -1(a |a |<1

14、向量 a 102-1z R =cot h -1(a |a |>1向量 a +1a -102z R =0. 5ln (a a >0向量 a +1/4a -1/40z R =ln (a/4 a >0向量a +ba -bz R =0. 5ln (a/b ab >03 C OR DIC 内核及前后处理单元设计由于在 CORDIC 迭代算法中如计算 32bit 具有 IEEE 2754标准格式的正余弦值需要 34个时钟周期 4, 其计算速 度很慢 , 对于大量的计算来说 , 不能够最大限度地利用 CPU 。 本设计采用流水线技术 5, 以面积换取速度的方式 , 通过增加硬件资源来

15、提高计算速度 , 使 CORDIC 的每一级 都同 时 工 作 。考 虑 到 实 际 计 算 复 杂 度 和 计 算 精 度 , CORDIC 内核中采用 32bit 宽的数据和 32级流水线结构 , 可在 1时钟周期计算出结果 , 延时 32个时钟周期而已 。由于迭代运算的次数有限 ,CORDIC 前置处理设计中 , 内核的数据表示为 :用 32bit 无符号二进制数表示浮点数 1。 对 0max 范围内的角度用 31bit 有符号二进制数来表 示 , 后置处理单元在 CORDIC 输出指数函数值时用 33bit 有符号二进制数整数来表示 , 以防数据溢出 。基于流水线技术的 CORDIC

16、硬件实现架构 2,3如图 2151 图 3 功能仿真波形 表 2 误差分析x 1/161/81/40. 3-1/8-0. 3x 二进制表示134923404-332344144相对误差2. 976e -8-1. 482e -8-0. 730e -8-1. 482e -8所示 。 图 2中架构由控制单元 、 选择器 、 移位寄存器 、 ROM和加法器组成 。 ROM 中 存 放 的 是 预 先 计 算好 的 32个tanh -1数据值 4:n =tanh -1(1/2n 3230, n =1, 2, , 32, 供处理器运算时调用 。 从架构中可以看出只有移位和加减法运算 , 因此特别适合于可编

17、程逻辑器件 。 每次计算之前 要对 Y 值进行判断 , 根据 Y 的符号来决定图 2中的加减运 算 4, 得到的结果再通过寄存器保存中间变量供下级流水 线使用 。由于将 CORDIC 运算由循环迭代展开为多级流 水线结构 , 级级直接相连 , 相当于直接连线 , 每级移位宽度 固定 , 减少了循环迭代中需要的桶型算术移位寄存器。 图 2 CORDIC 流水线处理器架构4 功能仿真及时序验证用 Modelsim 软件对本设计进行仿真 , 选取 1/16、 1/8、1/4、 0. 3、 -1/8、 -0. 3等 6个值进行功能仿真 , 如图 3所示 。 从仿真图看出输入值 iAngle 为 -0.

18、 3(十进制表示为-322122547 , 经过延时 32时钟周期计算出来 , 由此可 知利用流水线技术可以不断地计算出反双曲正切函数 , 其 它值依此类推 。 在 FP GA 实现反双曲正切函数的计算 , 对 于需要精确计算反双曲正切函数的应用场合来说显得非常 重要 。表 2给出了不同数值的功能仿真误差分析 ,tanh -1(1/16 的绝对误差 atanh =2, 为 atanh 计算出来的实际值与 理论值 之 差 , 对 应 的 相 对 误 差 为 atanh/274247418=2. 976e -8。 从表 2的误差分析可知 , x 输入在一定范围内反双曲正切函数的绝对误差很小 , 相

19、对误差的数量级为10-8, 其精度非常高 , 可以满足很多高精度计算场合的需要 。 例如 ,LDPC 码在对数似然比域和最小和算法 BP (置信传播算法 译码的 FP GA 实现 , 需要在译码过程中使用 反双曲正切函数来计算校验节点信息的更新 6, 其译码计 算过程要求具有实时性 、 高精度性以保证译码的准确性 、 可 靠性 。 反双曲正切函数对 LDPC 码的译码具有很重要的意 义 , 为 LDPC 译码的硬件实现奠定了基础 。在 Quartus II 6. 0软件上进行时序验证 , 其时序结果 与功能仿真的结果一致 , 并在 Altera 公司的 C YCLON E II EP2C356

20、72C6器件上综合 , 消耗 4437个 L E (逻辑单元 、 2911个 REG (寄存器 , 运行 Timing Analyzer Tool (时序分 析工具 , 最高的工作频率为 161. 45M Hz , 并在 DE2开发 板上进行硬件下载验证 , 该算法可以满足 实际的需求 , 。32位 反 双 曲 正 切 函 数 是 采 用 流 水 线 技 术 CORDIC 算法来实现的 , 以面积换速度 , 能够在 1个时钟周 期内计算出函数值 , 达到了很高的精度和很快的速度 , 完成 了反双曲正切函数的 IP 核部分 , 具有较广泛的应用价值 。 也可以在 Nios 处理器中增加自定义指令 , 即将一