前三篇文章里我用定点查表法把sinf()的计算速度提升了14倍把atan2f()提升了8倍把sqrtf()提升了2.46倍。今天我们来手撕expf()——这个在数字滤波、概率计算、衰减模型里无处不在的函数。01 为什么是exp在嵌入式开发中exp的出场频率比你想象的高- 数字滤波一阶低通滤波器的系数计算- 传感器线性化热敏电阻、光电二极管的非线性校正- 概率计算贝叶斯推断、softmax- 衰减模型电池放电曲线、RC充放电- 音频处理指数包络、动态范围压缩标准库的expf()在Cortex-M4上虽然可以用FPU加速部分运算但内部有大量分支和范围检查。然而实际工程应用中输出通常在某个固定范围之内根据实际情况制定定点fast_exp()用拆分整数小数加查表和插值的方式可以做到更快、更可控。02 exp的特殊难点动态范围极大和sin()、atan2()、sqrt()不同exp()有一个天然的难点它的值域跨度极大。输入 x-10-101510exp(x)0.0000450.3681.02.718148.422026Q15 格式只能表示-1.0 ~ 0.99997根本装不下exp()的全部结果。所以定点exp() 不能像sin()那样直接查表输出必须做范围压缩。03 核心思路拆整数小数指数函数有一个天然性质其中- k 是整数部分- r 是小数部分范围 |r|≤ln2/2 ≈ 0.3466- 2^k 用移位实现- exp(r) 在[-0.3466, 0.3466] 范围内值域约为[0.707, 1.414]Q15 装不下 1.414关键改进我们不用四舍五入取k而是用ceil 向上取整这样exp(r) in (0.5, 1.0]Q15 刚好装得下不需要任何额外的缩放。所以流程是输入 x (Q6.10)计算 k ceil(x / ln2)计算 r x - k·ln2 (范围 (-ln2, 0])查表求 exp(r) (Q15, 范围 (0.5, 1.0])输出尾数 exp(r) 和指数 k实际值 尾数 / 32768.0 × 2^k关键洞察exp()的定点实现不是查表求exp()而是查表求exp(r) 移位求2^k。移位是零成本的所以真正耗时的只有查表那一步。04 fast_exp实现exp_table[258]exp(r) 的 Q15 值r (i - 256) / 256i 0~256extern const uint16_t exp_table[258];// 常数Q15 格式#define INV_LN2_Q15 47274 // (1/ln2) × 2^15 ≈ 47274#define LN2_Q15 22713 // ln2 × 2^15 ≈ 22713typedef struct {int16_t mantissa; // Q15范围 [0.3679, 1.0]int16_t exponent; // 2 的幂次} exp_result_t;/*** brief 快速指数函数* param x 输入Q6.10 格式实际值 x / 1024* return exp_result_t实际值 mantissa / 32768.0 × 2^exponent** 推荐输入范围x ∈ [-32768, 32767]实际值 [-32,31.999]* /exp_result_t fast_exp(int16_t x) {// 1. 计算 t x / ln2Q6.10 格式乘法代替除法四舍五入int32_t t ((int32_t)x * INV_LN2_Q15 16384) 15;// 2. 计算 k ceil(t / 256)无分支向上取整int32_t k (t 1023) 10;// 3. 计算 r x - k × ln2Q6.10 格式int32_t r_q6_10 (int32_t)x - ((k * LN2_Q15 16) 5)1024;int32_t idx r_q6_10 2;int32_t frac r_q6_10 0x03;// 4. 查表exp_result_t result;int32_t res (int32_t)exp_table[idx1] - exp_table[idx];res (res * frac 2)2;res exp_table[idx];result.mantissa (int16_t)res;result.exponent (int16_t)k;return result;}关于除法代码中用乘法 x * INV_LN2_Q15 代替了除法。在 Cortex-M4 上乘法只需 1~2 个周期比除法快得多。如果芯片没有硬件除法如 M0这个优化更加明显。关于执行时间查表移位的方式没有循环每次调用的执行时间都是固定的对实时控制非常友好。05 使用示例// 计算 exp(1.0)// x 1.0 × 1024 1024exp_result_t r fast_exp(1024);// 结果r.mantissa ≈ 22259r.exponent 2// 实际值 22259 / 32768.0 × 2^2 ≈ 2.71716// 计算 exp(-0.5)// x -0.5 × 1024 -512exp_result_t r2 fast_exp(-512);// 结果r2.mantissa 19875r2.exponent 0// 实际值 19875/32768 ≈ 0.60654// 计算 exp(0)exp_result_t r3 fast_exp(0);// 结果r3.mantissa 32767r3.exponent 0// 实际值 32767 / 32768.0 × 1 ≈ 0.9999706 精度分析fast_exp 的误差主要来自三个环节误差来源量级说明查表插值0.04%256点查表r∈(-ln2,0]Q15尾数量化0.003%表格值存储为Q15输入Q6.10量化0.05%输入x的分辨率1/1024总相对误差0.05%-0.09%主要瓶颈是输入量化注意exp() 的误差应该用相对误差衡量因为输出值域跨度极大从1.24x10^-14到7.9x10^13绝对误差没有意义。如果要降低总误差需要提高输入格式的位宽例如改用 Q16.16。07 性能实测测试平台GD32F303VGT6 120MHzDWT周期计数32768次调用。函数开启FPU周期关闭FPU周期expf()366994720164341fast_exp()2228224222822408 源码与工程完整的Keil工程已上传Gitee包含- GD32F303VGT6的完整工程- fast_exp() 完整源码- 性能测试代码使用DWT周期计数- 精度测试代码Gitee地址https://gitee.com/jervis_luo/gd32_fast_exp.git09 总结对比维度expf()开FPUexpf()关FPUfast_exp()32768次耗时3,669,94720,164,3412,228,224最大相对误差~1e-7~1e-70.05%-0.09%是否需要FPU是否否适用芯片M4/M7M0/M3/M4/M7M0/M3/M4/M7三点核心收获1. 定点exp()的核心是范围压缩——exp(x) 2^k · exp(r)用ceil取整保证r≤0让exp(r)落在Q15范围内。2. 误差用相对误差衡量——exp的值域跨度极大绝对误差没有意义总误差主要来自输入量化。3. 执行时间固定——没有循环和分支预测失败对实时控制非常友好。系列文章1. [fast_sin比官方sinf快14倍]2. [fast_atan2FOC角度计算优化128点查表法让atan2f快8倍]3. [fast_sqrtFOC矢量幅值计算优化256点查表牛顿迭代让sqrtf快2.5倍]4. [手撕官方expf拆整数查表插值关闭FPU快9倍](本文)后记写完这篇文章后我把 fast_exp() 用在了自己做的电池管理模块里。原本用 expf() 计算SOC曲线时每个采样周期要花几十微秒换上 fast_exp() 后直接降到几微秒。整个模块的采样率从1kHz提升到了10kHz。这大概就是优化的乐趣所在。如果你对定点数学库感兴趣欢迎关注我。这个系列到这里已经覆盖了sin()、atan2()、sqrt()、exp() 四个函数后续如果大家感兴趣我可以继续写 fast_log()、fast_pow() 或者矩阵运算的定点实现。