从纸笔到芯片:手把手拆解CPU除法器的前世今生(附Verilog代码)
在计算机体系结构的浩瀚星河中,除法器始终是那颗既令人着迷又让人望而生畏的星辰。当我们用Python写下简单的a/b时,CPU内部究竟上演着怎样的微观戏剧?从古希腊人用算盘计算比值,到现代超标量处理器每个时钟周期完成多次除法,这段跨越两千年的技术进化史,隐藏着人类对计算效率永无止境的追求。
本文将带你穿越时空长廊,从最基础的纸笔算法出发,逐步拆解六种关键除法器设计范式。我们会用Verilog代码具象化每个阶段的突破,最终在FPGA上实现一个完整的基2 SRT除法单元。不同于教科书式的概念罗列,我们更关注为什么某种设计会在特定历史节点胜出——是晶体管预算的制约?还是时钟频率的瓶颈?或是并行计算的需求?
1. 计算之始:纸笔算法的硬件映射
当我们用竖式计算124÷3时,大脑其实在执行一套精密的决策流程:比较部分余数与除数、试探性相乘、结果修正。这套源自印度的算法(史称"galley division")在1947年被ENIAC的设计者们首次搬进电子计算机。
核心迭代公式:
R_{i+1} = B × R_i - q_{n-i-1} × D其中B表示基数(十进制为10,二进制为2),D是除数,q是当前位的商。这个看似简单的等式,将在后续所有算法变种中反复出现。
纸笔算法的Verilog实现暴露了根本缺陷:
module paper_pencil_divide #(parameter WIDTH=8) ( input [WIDTH-1:0] dividend, input [WIDTH-1:0] divisor, output reg [WIDTH-1:0] quotient, output reg [WIDTH-1:0] remainder ); integer i; reg [2*WIDTH-1:0] partial_remainder; always @(*) begin partial_remainder = dividend; quotient = 0; for(i=0; i<WIDTH; i=i+1) begin if(partial_remainder >= (divisor << (WIDTH-1-i))) begin quotient[WIDTH-1-i] = 1'b1; partial_remainder = partial_remainder - (divisor << (WIDTH-1-i)); end end remainder = partial_remainder; end endmodule这个实现需要WIDTH个周期完成计算,关键路径包含:
- 移位器(divisor << n)
- 比较器(partial_remainder ≥ ...)
- 减法器
在早期电子管计算机中,这种串行设计导致除法耗时是乘法的10倍以上。这催生了我们对下一个里程碑的探索——如何通过硬件并行化加速?
2. 机械革命:恢复余数法的硬件实现
1951年,冯·诺伊曼在IAS计算机中首次应用恢复余数法(Restoring Division),其核心创新在于用硬件状态机替代人工决策。算法流程如下:
- 初始化:将除数左移n位,被除数存入A寄存器
- 迭代步骤:
- 执行A = A - D
- 若结果为负:商位置0,执行A = A + D(恢复)
- 若结果为正:商位置1
- 左移A寄存器
Verilog实现揭示关键优化:
module restoring_divide #(parameter WIDTH=4) ( input clk, start, input [WIDTH-1:0] dividend, divisor, output reg [WIDTH-1:0] quotient, output reg [WIDTH-1:0] remainder, output reg busy ); reg [2*WIDTH:0] A; // 额外1位用于符号判断 reg [2*WIDTH:0] D; reg [WIDTH-1:0] counter; always @(posedge clk) begin if(start) begin A <= { {WIDTH{1'b0}}, dividend, 1'b0 }; D <= { divisor, {WIDTH+1{1'b0}} }; quotient <= 0; counter <= WIDTH; busy <= 1; end else if(busy) begin // 试探性减法 A <= A - D; if(A[2*WIDTH]) begin // 结果为负 quotient <= { quotient[WIDTH-2:0], 1'b0 }; A <= A + D; // 恢复 end else begin quotient <= { quotient[WIDTH-2:0], 1'b1 }; end D <= D >> 1; counter <= counter - 1; if(counter == 0) begin remainder <= A[2*WIDTH:WIDTH+1]; busy <= 0; end end end endmodule性能瓶颈分析:
| 操作 | 延迟(时钟周期) | 关键路径 |
|---|---|---|
| 减法器 | 1 | 全位宽并行减法 |
| 结果判断 | 1 | 最高位符号检测 |
| 恢复操作 | 1 | 多路选择器+加法器 |
| 总计/位 | 3 | 最长组合逻辑链 |
虽然比纸笔算法快了一个数量级,但每比特商需要3个时钟周期的设计在GHz时代显得过于奢侈。这引出了现代除法器的第一个重要范式转变——消除恢复操作。
3. 速度飞跃:不恢复余数法的精妙设计
1957年,IBM 704首次采用不恢复余数法(Non-Restoring Division),其革命性在于允许部分余数暂时为负,通过后续操作自动修正。算法步骤重构为:
- 若部分余数为正:商1,下轮迭代执行减法
- 若部分余数为负:商0,下轮迭代执行加法
Verilog实现展现面积优化:
module nonrestoring_divide #(parameter WIDTH=8) ( input clk, start, input [WIDTH-1:0] dividend, divisor, output reg [WIDTH-1:0] quotient, output reg [WIDTH-1:0] remainder, output reg busy ); reg [2*WIDTH:0] A; reg [2*WIDTH:0] D; reg [WIDTH-1:0] counter; reg last_negative; always @(posedge clk) begin if(start) begin A <= { {WIDTH{1'b0}}, dividend }; D <= { divisor, {WIDTH{1'b0}} }; quotient <= 0; counter <= WIDTH; busy <= 1; last_negative <= 0; end else if(busy) begin if(A[2*WIDTH]) begin // 当前余数为负 A <= { A[2*WIDTH-1:0], 1'b0 } + D; quotient <= { quotient[WIDTH-2:0], 1'b0 }; end else begin A <= { A[2*WIDTH-1:0], 1'b0 } - D; quotient <= { quotient[WIDTH-2:0], 1'b1 }; end counter <= counter - 1; last_negative <= A[2*WIDTH]; if(counter == 0) begin // 后处理:修正最终余数 if(last_negative) begin remainder <= A[2*WIDTH:WIDTH] + divisor; quotient <= quotient - 1; end else begin remainder <= A[2*WIDTH:WIDTH]; end busy <= 0; end end end endmodule关键改进对比:
| 特性 | 恢复余数法 | 不恢复余数法 |
|---|---|---|
| 迭代周期/bit | 3 | 1 |
| 硬件面积 | 加法器+减法器 | 单加法器+符号控制 |
| 关键路径延迟 | 减法+恢复 | 条件加法 |
| 功耗 | 高(冗余操作) | 低(最小化操作) |
这种设计在1960年代的CDC 6600超级计算机中表现惊艳,将除法速度提升至与乘法相当。但当RISC架构兴起时,人们发现新的瓶颈——商位选择逻辑过于串行化。
4. 现代基石:SRT算法的并行化突破
1989年,IBM研究员Sweeney、Robertson和Tocher提出的SRT算法(以三人姓氏首字母命名)带来三大革新:
- 冗余数字集:允许商位选择存在临时不确定性
- 提前终止比较:仅需检查部分余数的最高几位
- 并行进位保存:使用CSA(Carry-Save Adder)加速迭代
基2-SRT的商位选择规则:
| 部分余数范围 | 商位选择 |
|---|---|
| P ≥ 0.5 | +1 |
| -0.5 ≤ P < 0.5 | 0 |
| P < -0.5 | -1 |
Verilog实现展示现代设计范式:
module srt_radix2 #(parameter WIDTH=32) ( input clk, reset, input [WIDTH-1:0] dividend, input [WIDTH-1:0] divisor, output reg [WIDTH-1:0] quotient, output reg [WIDTH-1:0] remainder, output reg done ); // 规格化处理 wire [WIDTH-1:0] normalized_divisor; wire [WIDTH-1:0] normalized_dividend; wire [5:0] shift_amount; normalizer norm_div ( .in(divisor), .out(normalized_divisor), .shift(shift_amount) ); normalizer norm_dividend ( .in(dividend), .out(normalized_dividend) ); // 迭代状态机 reg [WIDTH+1:0] partial_remainder; reg [WIDTH-1:0] div_reg; reg [2*WIDTH-1:0] quotient_reg; reg [5:0] iteration; reg running; always @(posedge clk or posedge reset) begin if(reset) begin running <= 0; done <= 0; end else if(!running && !done) begin // 初始化 partial_remainder <= { {2{1'b0}}, normalized_dividend }; div_reg <= normalized_divisor; quotient_reg <= 0; iteration <= WIDTH; running <= 1; end else if(running) begin // 商位选择 casez(partial_remainder[WIDTH+1:WIDTH-2]) 4'b01??: begin // P ≥ 0.5 partial_remainder <= { partial_remainder[WIDTH-1:0], 2'b00 } - { div_reg, {WIDTH{1'b0}} }; quotient_reg <= { quotient_reg[2*WIDTH-3:0], 2'b01 }; end 4'b000?,4'b001?: begin // -0.5 ≤ P < 0.5 partial_remainder <= { partial_remainder[WIDTH-1:0], 2'b00 }; quotient_reg <= { quotient_reg[2*WIDTH-3:0], 2'b00 }; end 4'b11??,4'b10??: begin // P < -0.5 partial_remainder <= { partial_remainder[WIDTH-1:0], 2'b00 } + { div_reg, {WIDTH{1'b0}} }; quotient_reg <= { quotient_reg[2*WIDTH-3:0], 2'b11 }; // -1的补码表示 end endcase iteration <= iteration - 1; if(iteration == 0) begin // 后处理 if(partial_remainder[WIDTH+1]) begin // 余数为负 remainder <= partial_remainder[WIDTH+1:2] + div_reg; quotient <= quotient_reg[2*WIDTH-1:WIDTH] - 1; end else begin remainder <= partial_remainder[WIDTH+1:2]; quotient <= quotient_reg[2*WIDTH-1:WIDTH]; end done <= 1; running <= 0; end end end endmoduleSRT算法的三大加速技巧:
- 商位预测表(QDS):通过预计算的2-3位比较结果直接确定商位
// 简化的QDS查找表 module qds_lut ( input [3:0] partial_remainder_msb, input [3:0] divisor_msb, output reg [1:0] quotient_digit ); always @(*) begin case({partial_remainder_msb, divisor_msb}) // P ≥ 0.5 8'b01??_????: quotient_digit = 2'b01; // +1 // -0.5 ≤ P < 0.5 8'b000?_????, 8'b001?_????: quotient_digit = 2'b00; // 0 // P < -0.5 default: quotient_digit = 2'b11; // -1 endcase end endmodule- 进位保存加法器(CSA):将三个数相加的延迟从2个全加器降为1个
module csa #(parameter WIDTH=4) ( input [WIDTH-1:0] a, b, c, output [WIDTH-1:0] sum, carry ); assign sum = a ^ b ^ c; assign carry = { (a & b) | (a & c) | (b & c), 1'b0 }; endmodule- 实时商转换(On-the-Fly):避免最终冗余到常规表示的额外周期
module on_the_fly_converter ( input clk, input [1:0] current_q, // 当前商位(-1,0,+1) output reg [WIDTH-1:0] quotient ); reg [WIDTH-1:0] Q_plus; // Q+版本 reg [WIDTH-1:0] Q_minus; // Q-版本 always @(posedge clk) begin case(current_q) 2'b01: begin // +1 Q_plus <= { Q_plus[WIDTH-2:0], 1'b1 }; Q_minus <= { Q_minus[WIDTH-2:0], 1'b0 }; end 2'b00: begin // 0 Q_plus <= { Q_plus[WIDTH-2:0], 1'b0 }; Q_minus <= { Q_minus[WIDTH-2:0], 1'b1 }; end 2'b11: begin // -1 Q_plus <= { Q_minus[WIDTH-2:0], 1'b1 }; // Q- + ulp Q_minus <= { Q_minus[WIDTH-2:0], 1'b0 }; end endcase end // 最终选择Q_plus作为输出 assign quotient = Q_plus; endmodule5. 性能巅峰:基4与高基算法的工程实践
当SRT算法将基数从2提升到4时,每次迭代能产生2位商,理论上速度翻倍。但这也带来新的设计挑战:
- 商位选择更复杂(基4通常选择{-2,-1,0,1,2})
- 部分余数范围重叠区更关键
- 需要更精确的规格化处理
基4 SRT的PD图分析:
横轴为规格化除数d,纵轴为缩放后的部分余数P,阴影区域表示商位选择重叠区
现代处理器如Apple M1采用基16甚至更高基算法,配合以下优化技术:
- 商位选择预测:使用6-8输入查找表提前确定商位
- 多路并行CSA:同时计算多个可能路径
- 错误容忍设计:允许临时错误在后续迭代修正
module radix4_srt #(parameter WIDTH=64) ( input clk, reset, input [WIDTH-1:0] dividend, input [WIDTH-1:0] divisor, output [WIDTH-1:0] quotient, output [WIDTH-1:0] remainder, output done ); // 规格化处理 wire [WIDTH-1:0] normalized_divisor = divisor << leading_zeros(divisor); wire [WIDTH-1:0] normalized_dividend = dividend << leading_zeros(dividend); // 迭代核心 reg [WIDTH+3:0] partial_remainder; reg [WIDTH-1:0] div_reg; reg [2*WIDTH-1:0] quotient_reg; reg [5:0] iteration; // 基4 QDS模块 wire [2:0] q_selection; qds_radix4 qds ( .P(partial_remainder[WIDTH+3:WIDTH-4]), .D(div_reg[WIDTH-1:WIDTH-4]), .q(q_selection) ); always @(posedge clk) begin if(reset) begin iteration <= WIDTH/2; partial_remainder <= { 4'b0, normalized_dividend }; div_reg <= normalized_divisor; end else begin case(q_selection) 3'b000: // +2 partial_remainder <= (partial_remainder << 2) - (div_reg << 1); 3'b001: // +1 partial_remainder <= (partial_remainder << 2) - div_reg; 3'b010: // 0 partial_remainder <= partial_remainder << 2; 3'b011: // -1 partial_remainder <= (partial_remainder << 2) + div_reg; 3'b100: // -2 partial_remainder <= (partial_remainder << 2) + (div_reg << 1); endcase // 商位寄存器更新 quotient_reg <= { quotient_reg[2*WIDTH-3:0], q_selection }; iteration <= iteration - 1; end end // 后处理模块 post_processing pp ( .final_remainder(partial_remainder), .final_quotient(quotient_reg), .divisor(div_reg), .quotient(quotient), .remainder(remainder) ); assign done = (iteration == 0); endmodule6. 未来之路:融合乘法器的近似计算
在深度学习等新兴负载中,精确除法正被近似计算替代。两种前沿方向值得关注:
Newton-Raphson迭代法:将除法转化为乘法问题
x_{n+1} = x_n * (2 - d * x_n)需要3-4次迭代即可达到单精度浮点精度
Goldschmidt算法:同时收敛分子分母
k = 1/d的近似值 while(误差>阈值){ a = a * k b = b * k k = 2 - b }
这些算法在现代GPU如NVIDIA Tensor Core中广泛使用,配合低精度乘法器实现超高吞吐。