HLS / Chisel 实现CORDIC算法双曲系统

CORDIC(坐标旋转数字算法)是一种计算三角、双曲和其他数学函数的有效方法。它是一种数字算法，每次运算均产生一次结果输出。这使我们能够根据应用需求调整算法精度；增加运算迭代次数可以得到更精确的结果。

CORDIC是只使用加法、减法、移位和查找表实现的简单算法，这种算法在FPGA中实现效率高，在硬件算法实现中经常用到。

CORDIC算法双曲系统原理

点$(cosh\beta ,sinh\beta )$ 定义了右半直角双曲线$x^2-y^2=1$，所以有：

$$\{\begin{array}{l}{x}_{\mathrm{P}}=\cosh \beta \\ y_{\mathrm{P}}=\sinh \beta \end{array}$$

如下图我们假设$\beta$为自定义角度概念，然后旋转$\theta$角度到下一个点的这个过程我们做了一个旋转操作。（注意实际的角度非此大小，我们只是引入旋转概念）

同时有关于双曲三角函数的计算公式如下：

$$\{\begin{array}{l}\cosh (\theta \pm \beta )=\cosh \theta \cosh \beta \pm \sinh \theta \sinh \beta \\ \sinh (\theta \pm \beta )=\sinh \theta \cosh \beta \pm \cosh \theta \sinh \beta \end{array}$$

则可以令$\cosh \beta ={\mathbf{X}}_{\mathrm{i}},{\mathbf{Y}}_{\mathrm{i}}=\sinh \beta$ ，从中可以看出，将向量$(X_i，Y_i)$旋转$\theta$角，得到一个新的向量$(X_{i+1}，Y_{i+1})$ ,即上图的旋转，那么有：

$$\{\begin{array}{l}{X}_{i+1}=\cosh (\theta +\beta )=\cosh \theta \cosh \beta +\sinh \theta \sinh \beta =X_i\cosh \theta +Y_i\mathrm{shnh}\theta \\ Y_{i+1}=\sinh (\theta +\beta )=\sinh \theta \cosh \beta +\cosh \theta \sinh \beta =X_i\sinh \theta +Y_i\cosh \theta \end{array}$$

其中 $\theta$为旋转角度,化为矩阵表达形式, 平面旋转定义如下:

$$\left[\begin{array}{l}{\mathbf{X}}_{\mathrm{i}+1}\\ {\mathbf{Y}}_{\mathrm{i}+1}\end{array}\right]=\left[\begin{array}{ll}\cosh \theta & \sinh \theta \\ \sinh \theta & \cosh \theta \end{array}\right]\left[\begin{array}{l}{\mathbf{X}}_{\mathrm{i}}\\ {\mathbf{Y}}_{\mathrm{i}}\end{array}\right]$$

使用迭代的方法, 旋转角度可在多步完成, 每一步的旋转只完成其中的一小部分,多步之后将会完成一个平面旋转。旋转等式如上定义，旋转方法如下：

• 从一个起始位置开始旋转
• 每次旋转的一个先验固定的角度
• 每次旋转都朝向目标位置前进，不断调整方向逼近旋转目标

硬件优化

为了方便FPGA硬件计算，将如上的乘法转化为移位运算，则我们可以提取 $\cosh {\theta }_n$ 因子得到单步旋转等式为:

$$\begin{array}{rl}\left[\begin{array}{l}{\mathbf{X}}_{n+1}\\ {\mathbf{Y}}_{n+1}\end{array}\right]=& \left[\begin{array}{ll}\cosh {\theta }_n& \sinh {\theta }_n\\ \sinh {\theta }_n& \cosh {\theta }_n\end{array}\right]\left[\begin{array}{l}{\mathbf{X}}_n\\ {\mathbf{Y}}_n\end{array}\right]=\\ & \cosh {\theta }_n\left[\begin{array}{cc}1& \tanh {\theta }_n\\ \tanh {\theta }_n& 1\end{array}\right]\left[\begin{array}{l}{\mathbf{X}}_n\\ {\mathbf{Y}}_n\end{array}\right]\end{array}$$

然后令${\theta }_n={\tanh }^{-1}\left(2^{-n}\right)$，则$\tanh {\theta }_n=(2^{-n})$可以使得乘法转化为移位计算

由于${\tanh }^{-1}\left(2^{-0}\right)=\infty$, 其迭代应该从 n=1开始。为了区分旋转角度，加上方向因子${\sum }_{n=1}^{\infty }{\sigma }_n{\theta }_n=\theta ,{\sigma }_n=\pm 1$, 即所有迭代角之和必须等于旋转角度。

$$\left[\begin{array}{l}{\mathbf{X}}_{n+1}\\ {\mathbf{Y}}_{n+1}\end{array}\right]=\cosh {\theta }_n\left[\begin{array}{cc}1& S_n{2}^{-n}\\ S_n{2}^{-n}& 1\end{array}\right]\left[\begin{array}{l}{X}_n\\ Y_n\end{array}\right]$$

再通过 N 次迭代后：

$$\begin{array}{c}\left[\begin{array}{l}{\mathbf{X}}_{\mathrm{j}}\\ {\mathbf{Y}}_{\mathrm{j}}\end{array}\right]=\prod _{\mathrm{n}=1}^N\cosh {\theta }_{\mathrm{n}}\left[\begin{array}{cc}1& {\mathrm{S}}_{\mathrm{n}}{2}^{-\mathrm{n}}\\ {\mathrm{S}}_{\mathrm{n}}{2}^{-\mathrm{n}}& 1\end{array}\right]\left[\begin{array}{l}{\mathbf{X}}_{\mathrm{i}}\\ {\mathbf{Y}}_{\mathrm{i}}\end{array}\right]=\\ {\mathbf{K}}^{\ast }\prod _{\mathrm{n}=1}^N\left[\begin{array}{cc}1& {\mathrm{S}}_{\mathrm{n}}{2}^{-\mathrm{n}}\\ {\mathrm{S}}_{\mathrm{n}}{2}^{-\mathrm{n}}& 1\end{array}\right]\left[\begin{array}{l}{\mathbf{X}}_{\mathrm{i}}\\ {\mathbf{Y}}_{\mathrm{i}}\end{array}\right]\end{array}$$

将比例因子从迭代公式中提取出来, 定义 K 为增益因子, 具体增益K 取决于迭代次数， 对于所有的初始 向量和所有的旋转角度而言, K 是一个常数。通常把 K 称作聚焦常数,常数P为K的倒数，则有:

$$\begin{array}{c}K=1/P=\prod _{n=1}^N\cosh \left({\tanh }^{-1}\left(2^{-n}\right)\right)=\\ \prod _{n=1}^{\infty }\sqrt{1-2^{-2n}}\approx 0.8297816\end{array}$$

则每一次迭代我们可以只计算如下的公式，最后再乘上K,注意如下出现的z是对于角度的变化：

$$\{\begin{array}{l}{x}_{i+1}=x_i+{\sigma }_i\cdot y_i\cdot 2^{-i}\\ y_{i+1}=y_i+{\sigma }_i\cdot x_i\cdot 2^{-i}\\ z_{i+1}=z_i-{\sigma }_i\cdot {\tanh }^{-1}{2}^{-i}\end{array}$$

同时可以注意到如下的旋转的角度由于极限所以会限制旋转的大小只能在如下角度范围内。

$${\theta }_{\max }=\sum _{j=1}^{\infty }{\tanh }^{-1}\left(2^{-j}\right)\approx 1.0554693737354848$$

计算模式

双曲系统有旋转模式和向量模式两种计算模式。

旋转模式中,(x,y)在(1,0)点，Z为初始化需要旋转的角度（即与目标点的差距）, 当Z旋转变为 0时,则旋转完成，用来计算双曲正余弦函数，其推导公式如下：

$$\begin{gathered} \begin{array}{rl}& \{\begin{array}{l}{x}_n=K_n\left(x_1\cosh z_1+y_1\sinh z_1\right)\\ y_n=K_n\left(y_1\cosh z_1+x_1\sinh z_1\right)\\ z_n=0\\ K_n={\prod }_{i=1}^{n-1}\sqrt{1-2^{-2i}}\end{array}\\ & {\sigma }_i=\{\begin{array}{l}+1,z_i⩾0\\ -1,z_i<0\end{array}\end{array} \\ {z}_0\in [-\frac{\pi }{3},\frac{\pi }{3}] \end{gathered}$$

向量模式跟旋转模式差不多,其区别是旋转的方向取决于Y而不是 Y的符号, 即将输入向量 Y变为 0 , 计算相应的 ${\tanh }^{-1}$ (反双曲正切)值。几何解释相当于起点在目标向量（x,y）处，z记录角度的，旋转到x轴就可以累计得到z为向量与x轴的夹角；最终旋转到x轴时，x坐标即此向量的模，即长度。

$$\begin{gathered} \begin{array}{rl}& \{\begin{array}{l}{x}_n=K_n\sqrt{x_1^2-y_1^2}\\ y_n=0\\ z_n=z_1+{\tanh }^{-1}\left(y_1/x_1\right)\\ K_n={\prod }_{i=1}^{n-1}\sqrt{1-2^{-2i}}\end{array}\\ & {\sigma }_i=\{\begin{array}{l}+1,y_i⩽0\\ -1,y_i>0\end{array}\end{array} \\ -0.784<tanh(z_n)=\frac{y_1}{x_1}<0.784 \end{gathered}$$

可以看到如上都有输入输出$x_0,y_0,z_0$限制:

• 旋转模式中，由于${\theta }_{\max }={\sum }_{j=1}^{\infty }{\tanh }^{-1}\left(2^{-j}\right)\approx 1.055$可以知道必然其旋转不会超过$[-1.055,1.055]$ ,即本输入一定限制$z_0\in [-\frac{\pi }{3},\frac{\pi }{3}]$
• 向量模式中，注意$-0.784<tanh(z_n)=\frac{y_0}{x_0}<0.784$

HLS实践向量计算模式

#define THETA_TYPE float
#define COS_SIN_TYPE float
#define NUM_ITERATIONS 50

// The cordic_phase array holds the angle for the current rotation
// The data type of cordic & the iterations number of cordic
THETA_TYPE cordic_phase_h[] = {0.5493061443340548, 0.25541281188299536, 0.12565721414045303, 0.06258157147700301, 0.03126017849066699, 0.01562627175205221, 0.007812658951540421, 0.003906269868396826, 0.0019531274835325498, 0.000976562810441036, 0.0004882812888051128, 0.0002441406298506386, 0.00012207031310632982, 6.103515632579122e-05, 3.05175781344739e-05, 1.5258789063684237e-05, 7.62939453139803e-06, 3.8146972656435034e-06, 1.907348632814813e-06, 9.53674316406539e-07, 4.768371582031611e-07, 2.38418579101567e-07, 1.192092895507818e-07, 5.960464477539069e-08, 2.980232238769532e-08, 1.4901161193847656e-08, 7.450580596923828e-09, 3.725290298461914e-09, 1.862645149230957e-09, 9.313225746154785e-10, 4.656612873077393e-10, 2.3283064365386963e-10, 1.1641532182693481e-10, 5.820766091346741e-11, 2.9103830456733704e-11, 1.4551915228366852e-11, 7.275957614183426e-12, 3.637978807091713e-12, 1.8189894035458565e-12, 9.094947017729282e-13, 4.547473508864641e-13, 2.2737367544323206e-13, 1.1368683772161603e-13, 5.684341886080802e-14, 2.842170943040401e-14, 1.4210854715202004e-14, 7.105427357601002e-15, 3.552713678800501e-15, 1.7763568394002505e-15, 8.881784197001252e-16, 4.440892098500626e-16, 2.220446049250313e-16, 1.1102230246251565e-16, 5.551115123125783e-17, 2.7755575615628914e-17, 1.3877787807814457e-17, 6.938893903907228e-18, 3.469446951953614e-18, 1.734723475976807e-18, 8.673617379884035e-19, 4.336808689942018e-19, 2.168404344971009e-19, 1.0842021724855044e-19, 5.421010862427522e-20, 2.710505431213761e-20, 1.3552527156068805e-20, 6.776263578034403e-21, 3.3881317890172014e-21, 1.6940658945086007e-21, 8.470329472543003e-22, 4.235164736271502e-22, 2.117582368135751e-22, 1.0587911840678754e-22, 5.293955920339377e-23, 2.6469779601696886e-23, 1.3234889800848443e-23, 6.617444900424222e-24, 3.308722450212111e-24, 1.6543612251060553e-24, 8.271806125530277e-25, 4.1359030627651384e-25, 2.0679515313825692e-25, 1.0339757656912846e-25, 5.169878828456423e-26, 2.5849394142282115e-26, 1.2924697071141057e-26, 6.462348535570529e-27, 3.2311742677852644e-27, 1.6155871338926322e-27, 8.077935669463161e-28, 4.0389678347315804e-28, 2.0194839173657902e-28, 1.0097419586828951e-28, 5.048709793414476e-29, 2.524354896707238e-29, 1.262177448353619e-29, 6.310887241768095e-30, 3.1554436208840472e-30, 1.5777218104420236e-30};

void cordic_hv(COS_SIN_TYPE x, COS_SIN_TYPE y, COS_SIN_TYPE &x_out, THETA_TYPE &z_out)
{
    // Set the initial vector that we wil rorate
    // current_cos = I; current = Q
    THETA_TYPE theta = 0;
    
    // This loop iteratively rotates the initial vector to find the
    // sine and cosine vlaue corresponding to the input theta angle
    for(int j = 1;j <= NUM_ITERATIONS; j++){
    
        // Save ther current_cos, so that it can be used in the sine calculation
        COS_SIN_TYPE temp_x = x;
        
        // Perform ther roation
        if(y < 0){
            x = x + y >> j;
            y = y + temp_x >> j;
            theta = theta - cordic_phase_h[j-1];
        } else{
            x = x - y >> j;
            y = y - temp_x >> j;
            theta = theta + cordic_phase_h[j-1];
        }
    }
    
    // Set ther final value
    x_out = x;
    z_out = theta
}

Chisel实践向量计算模式

定点数定义

首先定义定点数，注意在chisel.experimental包里有定点数的类

import chisel3._
import chisel3.util._
import chisel3.experimental._
import scala.collection.immutable
import scala.math._

/* 定点数的位宽定义 */ 
trait HasDataConfig {
  val DataWidth = 24
  val BinaryPoint = 10
}

源码实现

然后实现cordic算法计算 输出x和arctanh：

• 其中注意有调试的接口将每一级流水线的计算值输出

import chisel3._
import chisel3.util._
import chisel3.experimental._
import scala.collection.immutable
import scala.math._


class cordic_hv extends Module with HasDataConfig with Cordic_Const_Data{
  /*
   * @x : 输入的向量横轴坐标 定点数类型
   * @y : 输入的向量横轴坐标 定点数类型 取值范围 y/x 属于(-1,1)
   * @x_out : 输出的x n次迭代的结果 定点数类型
   * @arctanh : 输出的arctanh(y/x)反三角函数 定点数类型
   * details: 利用cordic双曲坐标系的向量模式迭代得到双曲反三角函数的近似值，注意输入坐标
              在第一四象限，迭代次数不超过50，在25次时，K`值的变化已经超过了float的范围
  **/
  val io = IO(new Bundle {
    val x: FixedPoint = Input(FixedPoint(DataWidth.W, BinaryPoint.BP))
    val y: FixedPoint = Input(FixedPoint(DataWidth.W, BinaryPoint.BP))
    val x_out: FixedPoint = Output(FixedPoint(DataWidth.W, BinaryPoint.BP))
    val arctanh: FixedPoint = Output(FixedPoint(DataWidth.W, BinaryPoint.BP))
    /* Debug */
    //    val sin_o: Vec[FixedPoint] = Output(Vec(NUM_ITERATIONS, FixedPoint(DataWidth.W, BinaryPoint.BP)))
    //    val cos_o: Vec[FixedPoint] = Output(Vec(NUM_ITERATIONS, FixedPoint(DataWidth.W, BinaryPoint.BP)))
    //    val theta_o: Vec[FixedPoint] = Output(Vec(NUM_ITERATIONS, FixedPoint(DataWidth.W, BinaryPoint.BP)))
  })

  /* 旋转度数表 */
  val temp_seq: immutable.Seq[Double] = immutable.Seq(
    0.5493061443340548, 0.25541281188299536, 0.12565721414045303, 0.06258157147700301, 0.03126017849066699,
    0.01562627175205221, 0.007812658951540421, 0.003906269868396826, 0.0019531274835325498, 0.000976562810441036,
    0.0004882812888051128, 0.0002441406298506386, 0.00012207031310632982, 6.103515632579122e-05, 3.05175781344739e-05,
    1.5258789063684237e-05, 7.62939453139803e-06, 3.8146972656435034e-06, 1.907348632814813e-06, 9.53674316406539e-07,
    4.768371582031611e-07, 2.38418579101567e-07, 1.192092895507818e-07, 5.960464477539069e-08, 2.980232238769532e-08,
    1.4901161193847656e-08, 7.450580596923828e-09, 3.725290298461914e-09, 1.862645149230957e-09, 9.313225746154785e-10,
    4.656612873077393e-10, 2.3283064365386963e-10, 1.1641532182693481e-10, 5.820766091346741e-11, 2.9103830456733704e-11,
    1.4551915228366852e-11, 7.275957614183426e-12, 3.637978807091713e-12, 1.8189894035458565e-12, 9.094947017729282e-13,
    4.547473508864641e-13, 2.2737367544323206e-13, 1.1368683772161603e-13, 5.684341886080802e-14, 2.842170943040401e-14,
    1.4210854715202004e-14, 7.105427357601002e-15, 3.552713678800501e-15, 1.7763568394002505e-15, 8.881784197001252e-16,
    4.440892098500626e-16, 2.220446049250313e-16, 1.1102230246251565e-16, 5.551115123125783e-17, 2.7755575615628914e-17,
    1.3877787807814457e-17, 6.938893903907228e-18, 3.469446951953614e-18, 1.734723475976807e-18, 8.673617379884035e-19,
    4.336808689942018e-19, 2.168404344971009e-19, 1.0842021724855044e-19, 5.421010862427522e-20, 2.710505431213761e-20,
    1.3552527156068805e-20, 6.776263578034403e-21, 3.3881317890172014e-21, 1.6940658945086007e-21, 8.470329472543003e-22,
    4.235164736271502e-22, 2.117582368135751e-22, 1.0587911840678754e-22, 5.293955920339377e-23, 2.6469779601696886e-23,
    1.3234889800848443e-23, 6.617444900424222e-24, 3.308722450212111e-24, 1.6543612251060553e-24, 8.271806125530277e-25,
    4.1359030627651384e-25, 2.0679515313825692e-25, 1.0339757656912846e-25, 5.169878828456423e-26, 2.5849394142282115e-26,
    1.2924697071141057e-26, 6.462348535570529e-27, 3.2311742677852644e-27, 1.6155871338926322e-27, 8.077935669463161e-28,
    4.0389678347315804e-28, 2.0194839173657902e-28, 1.0097419586828951e-28, 5.048709793414476e-29, 2.524354896707238e-29,
    1.262177448353619e-29, 6.310887241768095e-30, 3.1554436208840472e-30, 1.5777218104420236e-30)
  val inits_cordic_phase_h: immutable.Seq[FixedPoint] = temp_seq.map(t => FixedPoint.fromDouble(t, DataWidth.W, BinaryPoint.BP))
  val cordic_phase_h: Vec[FixedPoint] = VecInit(inits_cordic_phase_h)

  /* 初始化计算的寄存器数组，形成NUM_ITERATIONS级流水 */
  val current_x: Vec[FixedPoint] = RegInit(VecInit(Seq.fill(NUM_ITERATIONS)(FixedPoint.fromDouble(0.0, DataWidth.W, BinaryPoint.BP)))) // cos
  val current_y: Vec[FixedPoint] = RegInit(VecInit(Seq.fill(NUM_ITERATIONS)(FixedPoint.fromDouble(0.0, DataWidth.W, BinaryPoint.BP)))) // sin
  val current_z: Vec[FixedPoint] = RegInit(VecInit(Seq.fill(NUM_ITERATIONS)(FixedPoint.fromDouble(0.0, DataWidth.W, BinaryPoint.BP)))) // 目标角度

  for (i <- 1 to NUM_ITERATIONS) {
    /*
     * x[i] = K`(x[i-1] + sigma * 2^(-i) * y[i-1])
     * y[i] = K`(y[i-1] + sigma * 2^(-i) * x[i-1])
     * z[i] = z[i-1] - sigma * arctanh(2^(-i))
    **/
    if (i == 1) {
      /* 流水线第一级直接对(io.x,io.y)点做运算 */
      when(io.y < FixedPoint.fromDouble(0.0, DataWidth.W, BinaryPoint.BP)) {
        current_x(i - 1) := io.x + (io.y >> i)
        current_y(i - 1) := io.y + (io.x >> i)
        current_z(i - 1) := -cordic_phase_h(i-1)
      }.otherwise {
        current_x(i - 1) := io.x - (io.y >> i)
        current_y(i - 1) := io.y - (io.x >> i)
        current_z(i - 1) := cordic_phase_h(i-1)
      }
    } else {
      when(current_y(i - 2) < FixedPoint.fromDouble(0.0, DataWidth.W, BinaryPoint.BP)) {
        current_x(i - 1) := current_x(i - 2) + (current_y(i - 2) >> i) // 移位替代乘法
        current_y(i - 1) := current_y(i - 2) + (current_x(i - 2) >> i)
        current_z(i - 1) := current_z(i - 2) - cordic_phase_h(i - 1)
      }.otherwise {
        current_x(i - 1) := current_x(i - 2) - (current_y(i - 2) >> i)
        current_y(i - 1) := current_y(i - 2) - (current_x(i - 2) >> i)
        current_z(i - 1) := current_z(i - 2) + cordic_phase_h(i-1)
      }
    }

    /* Debug */
    //    io.cos_o(i) := current_x(i)
    //    io.sin_o(i) := current_y(i)
    //    io.theta_o(i) := current_theta(i)
  }

  io.x_out := current_x(NUM_ITERATIONS - 1)
  io.arctanh := current_z(NUM_ITERATIONS - 1)
}

定义伴生对象工厂方法

最后我们定义两个个伴生对象工厂方法来方便直接调用函数，无需多加连线

/*
 * 定义这个类的伴生对象，并定义一个工厂方法来简化模块的例化和连线。
**/
object cordic_hv{
  def apply(x: FixedPoint, y: FixedPoint):(FixedPoint, FixedPoint) = {
    /*
    * @x : 输入的向量横轴坐标 定点数类型
    * @y : 输入的向量横轴坐标 定点数类型 取值范围 y/x 属于(-1,1)
    * @x_out : 输出的x n次迭代的结果 定点数类型
    * @arctanh : 输出的arctanh(y/x)反三角函数 定点数类型
    * details: 利用cordic双曲坐标系的向量模式迭代得到双曲反三角函数的近似值，注意输入坐标
    *          在第一四象限，迭代次数不超过50，在25次时，K`值的变化已经超过了float的范围
    **/
    val cordic_cv = Module(new cordic_hv)
    cordic_cv.io.x := x
    cordic_cv.io.y := y
    (cordic_cv.io.x_out, cordic_cv.io.arctanh)
  }
}
object cordic_arctanh{
  def apply(x: FixedPoint, y: FixedPoint):FixedPoint = {
    /*
    * @x : 输入的向量横轴坐标 定点数类型
    * @y : 输入的向量横轴坐标 定点数类型 取值范围 y/x 属于(-1,1)
    * @arctanh : 输出的arctanh(y/x)反三角函数 定点数类型
      details: 利用cordic双曲坐标系的向量模式迭代得到双曲反三角函数的近似值，注意输入坐标
               在第一四象限，迭代次数不超过50，在25次时，K`值的变化已经超过了float的范围
    **/
    val (x_out, arctanh) = cordic_hv(x, y)
    arctanh
  }
}

执行

object cordic_hv_app extends App {
  println(getVerilogString(new cordic_hv))
  (new chisel3.stage.ChiselStage).emitVerilog(new cordic_hv)
}

测试

import org.scalatest._
import chisel3._
import chisel3.experimental.BundleLiterals._
import chiseltest._
import chisel3.experimental._
import scala.math._

class Cordic_HV_Tester extends FlatSpec with HasDataConfig with ChiselScalatestTester with Matchers {
  behavior of "mytest2"
  it should "do something" in{
    val NUM_ITERATIONS = 20
    test(new cordic_hv){ c =>
      c.io.x.poke(FixedPoint.fromDouble(cosh(70*3.1415926/180),32.W, 20.BP))
      c.io.y.poke(FixedPoint.fromDouble(sinh(70*3.1415926/180),32.W, 20.BP))
      for(i <- 0 until NUM_ITERATIONS*2){
        println(s"*****************${i}******************")
        println(s"x_out ${c.io.x_out.peek}")
        println(s"arctanh ${c.io.arctanh.peek}\n")
        //         for(j <- 0 until NUM_ITERATIONS){
        //            println(s"sin_o  (${j}) ${c.io.sin_o(j).peek}")
        //            println(s"cos_o  (${j}) ${c.io.cos_o(j).peek}")
        //            println(s"theta_o(${j}) ${c.io.theta_o(j).peek}")
        //         }

        c.clock.step(1)
      }
    }
  }
}

参考论文

[1]陈石平, 李全, 莫丽兰,等. 基于CORDIC的反双曲正切函数的FPGA实现[J]. 计算机工程与科学, 2009, 31(5):3.