化石能源的获取、生产和分配等方面存在诸多问题，促使人们对可再生能源技术不断变革。随着电力电子技术和现代控制理论的发展，新能源发电使用的逆变器、新能源汽车使用的变频器技术迭代加速，这些领域的闭环实时控制技术也受到更多的关注，这对微控制器(microcontroller unit，MCU)的计算性能提出更高的要求。

变频器和逆变器的MCU控制算法的核心要素之一在于快速完成环路控制算法参数计算，从而保证控制系统的实时性和控制精度。PID控制的Park/Clark变换^[1]、NLPID算法中的输入整型^[2]、电机旋转变压器解码角度定位^[3]等常见的控制场景都需要计算三角函数、对数、指数函数等超越函数。另外，随着技术的迭代，较为先进的控制算法和更复杂的多电平拓扑结构都需要更高频次的超越函数计算能力^[4]。但是，使用通用CPU执行逼近算法对上述函数进行计算会消耗大量时钟周期，延长CPU的响应时间。在这种情况下，开发者往往会选择降低算法复杂度或采样频率从而满足实时性要求，这种方法虽能保持系统稳定性，却牺牲了控制精度。

超越函数的计算可以采用硬件加速器实现的方案降低运算周期^[5]，计算方法主要包含查表法、多项式拟合法和迭代法。单一查表法是通过映射输入数据和存储地址之间的关系，直接输出计算结果的一种方式，消耗芯片资源最少。但随着精度的提升和定义域的拓宽，表点的数量会呈现指数上升，因此只适用于低精度和小范围的计算^[6]。多项式拟合法通过构造特定形式的多项式，实现区间内的函数计算，主要包含泰勒级数展开、插值法等。多项式拟合法对于特定的函数区间具有很好的拟合精度，因此有些超越函数在较宽的输入区间会采用查表分段和拟合的方式进行组合处理^[7]，进而减少表点的数量。但是这种组合方式对于实现高精度的正弦/余弦计算的乘法器、加法器和查表单元来说，芯片资源消耗仍然相当可观^[8]。迭代法是硬件实现常用的一类方法，随着迭代次数的增加，计算精度具有线性收敛特性^[9]，常见算法包括牛顿拉弗森迭代算法、CORDIC(coordinate rotation digital computer)算法等。

作为数字迭代算法的一种，CORDIC算法可以以一种简洁的方式实现包括三角函数在内的多种超越函数的计算。该算法最早由Volder等^[10]于1959年提出，主要使用场景为美国航空控制系统，是一种基于平面直角坐标系旋转的算法。1971年Walther^[11]提出线性系统和双曲线系统统一的CORDIC算法，增加了乘除、双曲函数等计算的实现。虽然现代通用CPU的高主频和强性能使得CORDIC算法有了更多的替代方案，但至今该算法在控制、图像处理和通信领域的一些场景中仍然起到不可或缺的作用^[12]。为减少迭代的次数增强计算性能，基于角度重编码的CORDIC算法^[13]，设置多路比较器实现的并行计算可以进一步减少迭代计算的周期数^[14]。除减少迭代以外，CORDIC算法的硬件实现具有多种优化方式，可以使用定点数或者浮点数做迭代计算便于CPU进行数据读写^[15]，或者应用在加速非线性算法^[16]和FFT算法^[17]中。

然而，在工业数字实时控制场景应用，比如电机变频器和三相逆变器的控制算法中，高通量的数据流计算的需求较少，但实时性、配置灵活性和周期可预测性需求程度较高，而单一CORDIC算法提供的指数和对数计算更换底数方案较为繁琐。面向工控MCU中的实际应用需求，本文基于超越函数的迭代算法特点，设计了一种2通道8线程的超越函数单元(transcendental function unit，TFU)，支持定点小数q1.15和q1.31两种数据格式，可以实现正余弦、反正切、求模、求相角、对数和指数计算，同时具有快速响应、计算周期可预测、无中断快速读取的特性。

1 算法原理和场景需求 1.1 经典CORDIC算法计算三角函数

CORDIC的原理是通过多次迭代计算完成二维平面坐标轴向量的旋转，从而逼近待求解的目标角度。相比于其他类型的算法，CORDIC仅通过移位和加减法即可实现三角函数的运算。经典的CORDIC包括旋转模式和向量模式2种计算方式^[10]。

旋转模式将位于X轴正半轴上的单位向量进行逐步递减且收敛的旋转，当单位向量的旋转角等于平面直角坐标系指定角度θ时，其坐标值对应cosθ和sinθ。

在平面直角坐标系中点(x₁, y₁)和原点构成的向量旋转到(x₂, y₂)，当旋转角度为θ的时候，满足

$ \left\{\begin{array}{l} x_2=x_1 \cos \theta-y_1 \sin \theta , \\ y_2=x_1 \sin \theta+y_1 \cos \theta . \end{array}\right. $

(1)

CORDIC算法定义每一次迭代计算的角度正切值收敛，则迭代表达式可写作

$ \left\{\begin{array}{l} x_{2}=\cos \theta\left(x_{1}-y_{1} \tan \theta\right)=\cos \theta\left(x_{1}-d_{1} y_{1} 2^{-1}\right), \\ y_{2}=\cos \theta\left(y_{1}+x_{1} \tan \theta\right)=\cos \theta\left(y_{1}+d_{1} x_{1} 2^{-1}\right) . \end{array}\right. $

(2)

式中：d表示旋转角的方向，逆时针旋转时其值为1。对于每次迭代计算而言，不计算模长的伪旋转表达式为

$ \left\{\begin{array}{l} x_{i}=x_{i-1}-d_{i} y_{i-1} 2^{-i}, \\ y_{i}=y_{i-1}+d_{i} x_{i-1} 2^{-i}, \\ z_{i}=z_{i-1}-d_{i} \tan ^{-1} 2^{-i} . \end{array}\right. $

(3)

进行m次迭代后，可以得到

$ \left\{\begin{array}{l} x^{(m)}=K\left(x^{(0)} \cos z-y^{(0)} \sin z\right), \\ y^{(m)}=K\left(y^{(0)} \cos z+x^{(0)} \sin z\right), \\ z^{(m)}=0 . \end{array}\right. $

(4)

其中：K称为模长矫正因子或伸缩系数，随着迭代次数的增加收敛于一个固定的值。此时令x的初值为1/K，y的初值为0，则迭代计算m次后可以同时获得z的正弦和余弦。

向量模式将平面直角坐标系中的某个指定向量向X轴旋转，最终该向量到达X轴。向量模式的迭代表达式和旋转模式一致，仅符号位的判定有所区别。经过m次迭代后，方程满足

$ \left\{\begin{array}{l} x^{(m)}=K\left(x^{2}-y^{2}\right)^{1 / 2}, \\ y^{(m)}=0, \\ z^{(m)}=z+\tan ^{-1}(y / x) . \end{array}\right. $

(5)

CORDIC算法的收敛域较小，由等式(1)可知

$ \begin{gather*} z \in\left[-\sum\limits_{i=1}^{m} \tan ^{-1} 2^{-i}, \sum\limits_{i=1}^{m} \tan ^{-1} 2^{-i}\right] \approx \\ {\left[-99.7^{\circ}, 99.7^{\circ}\right], } \end{gather*} $

(6)

则对于平面直角坐标系中的任意角度，由于sin/cos为周期函数，因此可以将角度映射到[0, π/2]区间。对于非周期的反正切函数而言，可以通过查找表的缩放实现定义域的扩展。

1.2 数字迭代法计算对数和指数

下面的方程定义了对数函数的迭代运算^[5]

$ \left\{\begin{array}{l} x^{(i+1)}=x^{(i)} c^{(i)}=x^{(i)}\left(1+d_{i} 2^{-i}\right), \quad d_{i} \in\{-1, 0, 1\}, \\ y^{(i+1)}=y^{(i)}-\log _{2} c^{(i)}=y^{(i)}-\log _{2}\left(1+d_{i} 2^{-i}\right) . \end{array}\right. $

(7)

其中$1+d_{i} 2^{-i}$ 可以从查找表中读出。通过d_i的选择让x逼近1，在第m次迭代运算可以得到

$ \left\{\begin{align*} x^{(m)} & =x^{(0)} \prod c^{(i)} \Rightarrow \prod c^{(i)} \approx 1 / x, \\ y^{(m)} & =y^{(0)}-\sum \log _{2} c^{i}=y^{(0)}-\log _{2} \prod c^{(i)} \\ & \approx y^{(0)}+\log _{2} x . \end{align*}\right. $

(8)

通过设置y的初值为0可以计算得到log₂x。从式(8)易知x的收敛域为

$ x \in\left[\frac{1}{\prod\left(1+2^{-i}\right)}, \frac{1}{\prod\left(1-2^{-i}\right)}\right] \approx[0.41, 3.45] 。$

类似地，定义迭代运算

$ \left\{\begin{array}{l}x^{(i+1)}=x^{(i)}-\log _{2} c^{(i)}=x^{(i)}-\log _{2}\left(1+d_{i} 2^{-i}\right), \\ y^{(i+1)}=y^{(i)} c^{(i)}=y^{(i)}\left(1+d_{i} 2^{-i}\right) \quad d_{i} \in\{-1, 0, 1\} .\end{array}\right. $

(9)

选择d_i的值使x趋近于0，经过m次迭代可以得到

$ \left\{\begin{array}{l} x^{(m)}=x^{(0)}-\sum \log _{2} c^{i} \approx 0, \\ y^{(m)}=y^{(0)} \prod c^{(i)}=y 2^{\log _{2} \Pi c^{(i)}} \approx y 2^{x} . \end{array}\right. $

(10)

上述迭代方程中x的收敛域为

$ \begin{align*} x & \in\left[\sum \log _{2}\left(1-2^{-i}\right), \left[\sum \log _{2}\left(1+2^{-i}\right)\right]\right. \\ & \approx[-1.24, 1.56] . \end{align*} $

(11)

1.3 工业控制场景需求

闭环PWM控制是当今工业控制系统中常见的技术，主要应用于电机驱动、新能源发电并网等场景。一种典型的单环PID电机控制算法流程如图 1所示。测试数据显示，采用FPU计算单个三角函数需要约90个时钟周期，在200 MHz主频下时间约0.45 μs。由于高精度工控ADC采样率一般为4 M，单次采样和转换耗时0.25 μs，因此三角函数的计算时间长度远大于ADC采样和转换的时间，导致CPU需要较长时间等待三角函数计算结果，压缩了环内其他算法的计算时间。

工业控制场景对于超越函数计算单元的需求主要包含如下几个方面：

1) 更短的计算时间：总线读写时间+TFU计算时间小于CPU软件计算的时间；

2) 可配置的精度选择：工控ADC的精度为12 bit或者16 bit，因此TFU在多数应用下仅需计算半精度浮点或者q1.15定点小数；

3) 单核多线程的独立访问：高性能工控软件线程间计算相互独立，TFU的多线程支持可以减少线程间记录状态带来的时间损耗；

4) 双核系统的独立访问：高性能工控通常采用双核架构便于分配并行任务，因此TFU需要满足并行任务的计算。

基于上述分析，选用数字迭代算法实现超越函数的运算。

2 TFU设计 2.1 TFU算子的实现

算子即算法中的最小运算单元，定义实现每一次迭代计算的模块为一个算子。TFU以全组合逻辑的形式实现CORDIC算子和指数对数(LOGEXP)算子。

如图 2(a)所示，一个CORDIC算子可以实现4阶CORDIC迭代计算，每阶迭代需要从外部输入一个逐步递减的角度数据，并根据上一阶的计算结果对本阶计算的正负性进行选择。指数对数算子的电路结构如图 2(b)所示。和CORDIC算子类似，每个算子可以实现4阶指数或者对数的运算。

如表 1所示，通过配置不同的输入参数，相同的CORDIC算子可以实现5种函数的计算，LOGEXP算子实现2为底的对数和指数运算。其中反正切函数ATAN和求相角函数ATAN2实质上是同一种运算，通过配置第2参数Y为1可以实现ATAN运算，但此时第1输出没有实际意义。

2.2 TFU加速器架构

工控领域常用的ARM-M系列处理器具有较强的多线程运算性能，为减少线程之间共用外设和总线资源产生的运算能力下降，TFU支持双通道计算，每个通道最高支持4线程数据接收，所有通道共享AHB slave接口。最终实现的TFU结构框图如图 3所示。

其中CSR表示配置寄存器(control/status register)，主要实现函数不同参数的灵活配置。WDATA和RDATA分别对应函数的输入和输出，均由深度为2级的FIFO实现。每个通道中，4个线程的数据通过AHB总线写入后，由ARB模块进行仲裁，决定运算优先级，仲裁的方式可以配置为round-robin或者固定优先级顺序。被仲裁命中的数据将会传递到FSM模块中进行数据迭代，由于每次迭代的运算电路相同，因此多次迭代可以公用一组CORDIC或者LOGEXP算子，分别计算三角函数和指数、对数。此外，接口寄存器还包含LOGEXP_LUT域，用于LOGEXP算子更换底数。

FSM控制的数据流如图 4所示。当AHB slave接口接收数据后，无需额外配置使能寄存器，如果识别到数据已完成传输，状态机会直接跳转到工作模式。数据首先会转移到数据过滤器中，将多余或者出现格式错误的数据清除，随后经过数据映射模块将数据映射到算法的收敛区间。每次迭代计算的数据结果会被保存到专用的寄存器中，最后一次迭代产生的数据会输入到响应寄存器中等待总线读取。

2.3 接口与时序

TFU顶层接口采用AMBA3.0 AHB slave接口，通过时序逻辑判断数据地址，将接口时序转换为寄存器读写时序连接到4个线程的寄存器上。但是，CPU需要通过总线轮询或者进入中断获取各个外设状态，这2种类型的状态获取均需要消耗大量时钟周期，延长CPU获取运算结果的时间。

借助AHB总线的Wait-phase特性可以解决上述问题。当总线访问对某一线程的RDATA域时，若TFU该线程没有完成计算，则拉低HREADY接口，使总线停止执行后续的命令。接口时序图如图 5所示，总线对RDATA(RD)地址的访问，会将WDATA地址的写入命令延迟，直至RDATA中的数据A被读取至HRDATA接口后继续执行。

除接口寄存器外，TFU还支持DMA和CPU中断触发接口，通过配置寄存器中的指定域对应功能。

2.4 寄存器设计

TFU每个线程的寄存器包含CSR、WDATA、RDATA和LOGEXP_LUT 4种，根据功能分为只读和读写2种类型。其中CSR寄存器位宽为32 bit，占用一个总线地址，其详细描述如表 2所示。由于TFU算子运算周期数较少，为进一步减少数据读写占用的时间，TFU采用可配置的灵活设计。配置CSR prec域中的参数，可以控制函数的迭代次数，对于精度要求较低的场景可以通过牺牲冗余的精度实现更快的计算。配置CSR res_num域，可以实现只读取双输出函数的主输出，完成运算流程，释放下一个函数的运算。

如表 3所示，WDATA和RDATA分别具有2组32 bit位宽的寄存器，分别存储2个输入和输出数据。数据格式需要满足定点小数q1.15或q1.31，即1 bit符号位，后跟15 bit或者31 bit尾数位。对于LOGEXP_LUT域，根据算法原理，TFU需要48个查找表值，分别支持数字迭代计算中的加或者减运算，其格式满足q1.31。在默认状态下，该底数为2。

2.5 函数和数据映射

TFU支持7种超越函数的运算，每种函数会读取最多2组输入数据(argument, ARG)寄存器的数据，计算结果导入到2组输出数据(response, RES)寄存器中。表 4列出了各个函数的输入和输出寄存器配置要求。

由于CORDIC和LOGEXP算子具有收敛域，对于域外数据或者不便处理的域内数据，函数计算需要通过数据映射的策略进行预处理。其中，ATAN函数需要对CSR寄存器中的Scale域进行配置，从而实现对定义域外的数据进行计算。表 5统计了TFU各个函数的算法收敛域、TFU映射域和TFU最大可以支持的定义域。其中收敛域是指迭代计算算法原理可以支持的数据范围，三角函数用角度表示；映射域是TFU算子计算时的数据范围；定义域是TFU数据寄存器支持读取的数据范围。需要注意的是，LOGEXP算子的定义域是对收敛域取子集，因此不需要映射域进行域扩展。

3 结果与数据统计 3.1 TFU计算精度

TFU可以通过配置CSR寄存器实现1~6时钟周期的计算，每个时钟周期完成4次迭代计算。图 6展示TFU在不同迭代次数下最大误差的统计图。随着迭代次数的增加，函数计算的精度随之提升，并具有较好的线性度，与算法原理一致。值得注意的是，由于求模运算的结果需要乘以精度为31 bit的系数K，因此迭代次数小于16次时精度显著高于其他6个函数。各个函数的最大误差如表 6所示。图 7展示在工控系统中常用的q1.15数据格式下，迭代周期为4时的计算误差统计。

3.2 TFU综合

TFU通过Synopsys Design Compiler工具进行综合。由于TFU只有一个时钟域，所以设定所有的输入和输出端口的最大延迟为半个时钟周期。此外，设定setup和hold的uncertainty参数为0.1 ns，transaction ratio为0.8。此次综合选取SMIC 40 nm eFlash平台的40 nm沟道超高速常规阈值电压(40 nm-channel very-high-speed regular threshold voltage, VHSC40-RVT)标准单元库，实现了200 MHz的时钟频率，实际占用面积301 074 μm²。以NAND2标准单元估算逻辑门数量约110 688。

为便于评估模块的规模，对TFU使用Xilinx Vivado工具进行FPGA综合。结果显示，在使用Virtex7系列FPGA时，TFU共占用LUT数量为18 021，寄存器数量为6 244。

3.3 性能比较

TFU是一种面向工控MCU的专用外设，具有三角函数和指数对数的运算能力。和国内外相似文献相比，TFU具有高并行度、短计算周期、灵活配置的特性，并且提供了可以更换底数的指数对数算子。对比内容如表 7所示。

4 结论

本文面向工控MCU的算法实时性要求，设计了一种超越函数计算外设TFU。借助CORDIC算法和数字迭代计算，TFU可以实现正余弦、反正切、求模长以及对数和指数计算，具有4~24 bit可配置的计算精度。此外，TFU支持q1.15和q1.31定点小数格式，具有计算精度可配置的特性，借助总线接口可实现无中断触发的结果读取。面向双核系统架构，TFU具有2通道8线程，和可配置的线程优先级，对数指数计算可换底数，满足多场景需求。TFU使用SMIC40 eFlash标准单元库综合可以达到200 MHz的工作频率，综合面积为301 074 μm²。