高速高精度低资源占用单精度浮点算法verilog实现

最新推荐文章于 2025-06-10 15:55:32 发布

原创最新推荐文章于 2025-06-10 15:55:32 发布 · 7.6k 阅读

62 ·

CC 4.0 BY-SA版权

文章标签：

#verilog #ALTERA #fpga #浮点数 #单精度

FPGA 专栏收录该内容

7 篇文章

订阅专栏

本文介绍了一种针对FPGA平台的单精度浮点乘法算法优化方案，通过对IEEE754标准下的浮点数进行巧妙的拆分与重组，显著减少了所需硬件资源的同时保持了良好的计算精度与速度。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

高速高精度低资源占用单精度浮点算法verilog实现

最近一个需要考虑在FPGA里面实现一个对图像进行处理的算法，有几十次迭代计算，并且每次迭代还有很多浮点和矩阵乘法运算。FPGA板子中的乘法器数量是有限的，估计了一下远远不够使用，于是打算不用官方的浮点乘除法运算自己编写。直接使用一个ALTERA altfp_mult IP核需要使用7个9*9（4个18*18）乘法器模块，并且fmax直接编译出来只有160-180左右，经过重写的一个乘法器模块只用4个9*9(2个18*18)乘法器模块，并且精度、速度都几乎没有降低，fmax直接编译出来也有两百多。（芯片型号：使用DE2-115开发板）。

也希望本博客对读者有所帮助，文中的错误还望大家指正，不胜感激，如果转载请说明出处。

1，单精度浮点格式

双精度double型需要使用8个字节，占用资源太多，暂不考虑在FPGA上的实现，主要考虑单精度float型的实现。单精度浮点数的IEEE 745标准如图1所示，即符号位s占1bit，指数位e占8bits，尾数m占23bits。表示的小数为(-1)^s*2^(e-127)*1.m，注意e是加了127的偏移同时尾数是1.x中的x。

图 1 IEEE 745 单精度浮点数

2，单精度浮点数乘法

两个单精度浮点数相乘如(-1)^s1*2^(e1-127)*1.m1 x (-1)^s2*2^(e2-127)*1.m2结果为:
(-1)^s1+s2*2^e1+e2-127*1.m1*1.m2。因此只需要对1.m1*1.m2进行乘法运算。在FPGA中实现乘法显然可以对1.m1*1.m2进行直接相乘，是24bits*24bits，需要4个18*18乘法器模块。而FPGA中乘法器是以9*9,18*18的方式进行组织，直接算24bits*24bits会产生较大的延迟。可考虑将x,y拆分成小于18的两段相乘再相加加快速度，但占用资源数不会变。

图 2 乘法算法分割

为了降低资源使用量，考虑对算式进行精简、拆分。后面两个乘数均表示为1.x和1.y。首先拆分如下：

1.x*1.y=(1+x)*(1+y)=1+x+y+x*y; （1）

考虑到x,y是23位数，因此对X进行一次分割，分割为x1(18bits)和x2(5bits)，（1）式整理为：

1.x*1.y = (1+x1+x2)*(1+y1+y2)

= 1+x1+x2+y1+y2+x1*y1+x1*y2+x2*y1+x2*y2

= 1+x+y+ x1*y1+x1*y2+x2*y1+x2*y2 (2)

注意x2，y2小于2^-18，两数相乘则小于2^-36,超出单精度浮点数位数表示范围，因此可以将其省略。（2）式化简为（3）：

1.x*1.y = 1+x+y+ x1*y1+x1*y2+x2*y1 (3)

这里需要3个18*18乘法器模块。同样由于精度可以对x1*y2+x2*y1再进行精简近似。将x1拆分成x3(9 bits)和x4(5 bits),对y1进行同样的操作。则

x1*y2=(x3+x4)*y2 = x3*y2+x4*y2 (4),

由于x4小于2^-9,y2小于2^-18,则x4*y2小于2^-27，同样小于单精度浮点数位数表示范围，将其省略，对y进行同样的操作。因此(4)化简为

x1*y2 = x3*y2 (5)

综上所述，可以将1.x*1.y表示为

1.x*1.y = 1+x+y+ x1*y1+ x3*y2+x2*y3 (6)

x1*y1为18bits*18bits，x3*y2和x2*y3为9bits*5bits需要2两个9*9乘法器模块,共使用4个9*9乘法器模块（2个18*18乘法器模块），比官方的IP核使用少一半的资源，同时性能几乎没有降低。同样，调用此乘法器的除法器模块资源占用量也会大大减少。

3，模块的verilog pipeline方式实现

Pipeline方法充分利用了FPGA的并行特性。在实现两个数相乘的过程的模块中，首先将两个数加载进来，第一个时钟同时实现x+y，x1*y1，x3*y2，x2*y3，其后的时钟的时钟则分别将这些结果相加。
模块代码如下（尚未添加对0的异常处理等），在测试代码testbench中添加了自己编写的这个模块和ALTERA 的ALTFP_MULT IP核并在modelsim 10.a中仿真，以检查两个结果的差异。测试的浮点数使用google搜索到的 IEEE 754 Converter （ http://www.h-schmidt.net/FloatConverter/IEEE754.html ）产生。
模块代码：

module fpmul(
	input [31:0] mdat1,
	input [31:0] mdat2,
	output [31:0] odat,
	input clk,
	input rst_n
);

reg s1,s2;
reg [7:0] e1,e2;
reg [8:0] x3,x4,y3,y4;
reg [4:0] x2,y2;

always @(posedge clk or negedge rst_n) begin
	if(~rst_n) begin
		{s1,e1,x3,x4,x2}<=32'd0;
		{s2,e2,y3,y4,y2}<=32'd0;
	end else begin
		{s1,e1,x3,x4,x2}<=mdat1; //加载操作数
		{s2,e2,y3,y4,y2}<=mdat2; //加载操作数	
	end
end
// clk1
reg s3;
reg s3_b1,s3_b2,s3_b3,s3_b4;//pipeline保存数据

reg [8:0] adde;
reg [7:0] adde_b1,adde_b2,adde_b3,adde_b4;//pipeline保存数据

reg [35:0] mul1; //几个乘法计算
reg [13:0] mul2,mul3;
reg [23:0] addxy1;
//clk2
reg [22:0] mul1_b1;
reg [14:0] addm1;
reg [33:0] addm2;
//clk3
reg [33:0] addm3;
//clk4
reg [33:0] addm4;
//clk5
reg [22:0] addm5;
//output signal
assign odat={s3_b4,adde_b4,addm5};

always @(posedge clk or negedge rst_n) begin
	if(~rst_n) begin
		//
	end else begin
		//clk1
		s3<=s1^s2;
		adde<=e1+e2;
		mul1<={x3,x4}*{y3,y4};
		mul2<=x3*y2;
		mul3<=y3*x2;
		addxy1<={x3,x4,x2}+{y3,y4,y2};

		//clk2
		s3_b1<=s3;
		adde_b1<=adde-9'd127;
		addm1<=mul2+mul3;
		addm2<=mul1[35:4]+{addxy1,9'b0}; //23bits,对齐 //后面截止部分可能导致损失,因此增加了多余的9位进行相加
		//clk3
		s3_b2<=s3_b1;	
		adde_b2<=adde_b1;
		addm3<=addm2+addm1[14:0];//对齐
		//clk4
		s3_b3<=s3_b2;
		adde_b3<=adde_b2;
		addm4<=addm3+{1'b1,23'd0,9'b0};
		//clk5
		s3_b4<=s3_b3;
		if(addm4[33])begin //根据结果判定是否移位:1.x*1.y结果可能大于1,但是小于4，需要右移
			addm5<=addm4[32:10];
			adde_b4<=adde_b3+8'd1;
		end else begin
			addm5<=addm4[31:9];
			adde_b4<=adde_b3;
		end
	end
end

endmodule

测试模块testbench代码：

module FP_Test(
);

reg clk,rst_n;
reg [31:0] mdat1,mdat2;
wire [31:0] odat;
wire [31:0] odat2;

initial begin
	clk<=1'b1;
	rst_n<=1'b1;
	mdat1<=0;
	mdat2<=0;
	#2;
	rst_n<=1'b0;
	#2;
	rst_n<=1'b1;
	#4;
	mdat1<=32'h3fc30f28; //1.5239
	mdat2<=32'h3fc30f28;
	#2;
	mdat1<=32'h40a7a9fc; //5.2395
	mdat2<=32'h40a7a9fc;
	#2;
	mdat1<=32'h4251954d; //52.3958
	mdat2<=32'h4251954d;
	#2;
	mdat1<=32'h4402fd52; //523.9581
	mdat2<=32'h4402fd52;
	#2;
	mdat1<=32'h4251954d; //52.3958
	mdat2<=32'h4402fd52; //523.9581
	#2;
	mdat1<=32'h3f07929f; //0.529581
	mdat2<=32'h3f07929f;
	#2;
	mdat1<=32'hbf07929f; //-0.529581
	mdat2<=32'h4402fd52; //523.9581
	#2;
	mdat1<=32'h4380650b; //256.7894
	mdat2<=32'h4380650b; 
	#2;
	mdat1<=32'h3da1ab4b; //0.07894
	mdat2<=32'h3da1ab4b; 
	#2;
	mdat1<=32'h3da1ab4b; //0.07894
	mdat2<=32'h3f07929f; //0.529581
	#2;
	mdat1<=32'h3da1ab4b; //0.07894
	mdat2<=32'h0; //0.529581
	#2;
	mdat1<=32'h3bbab9a5; //0.0056984
	mdat2<=32'h3bbab9a5; //0.0056984
	#2;
	mdat1<=32'h3bbab9a5; //0.0056984
	mdat2<=32'h3da1ab4b; //0.07894

	#4;
	
end
always begin 
	#2;
	$display("sf:%x,alt:%x",odat,odat2);
end

always begin 
	#1;
	clk<=~clk;
end


fpmul mul0(
	.mdat1(mdat1),
	.mdat2(mdat2),
	.odat(odat),
	.clk(clk),
	.rst_n(rst_n)

);
alt_fpmul mul1(
	.aclr(~rst_n),
	.clock(clk),
	.dataa(mdat1),
	.datab(mdat2),
	.result(odat2),
	.zero());

endmodule

4，性能分析

在testbench中给了一些结果进行运算并仿真，这两个模块均只需间隔6个时钟即可得到结果。仿真波形如图3所示。

图 3 modelsim仿真波形

计算结果如下：

仅仅使用了这几个例子进行测试，可以发现两者结果基本一致，结果中我的模块有时候会比ALTERA IP核少1个bit的差异，会造成轻微的精度损失，这个损失基本上是可以接受的。需要注意的对于IEEE 745标准来说，整数部分越大，有效小数部分越小。
资源使用实际编译结果：

图4 使用ALTERA IP

图5 使用自己的模块

编译使用的顶层代码：

module fp_test1(
	input [31:0] mdat1,
	input [31:0] mdat2,
	output [31:0] odat1,
	input clk,
	input rst_n
);
reg [31:0] dat1,dat2,dat3;
always @(posedge clk or negedge rst_n) begin
	if(~rst_n) begin
		dat1<=0;
		dat2<=0;
		dat3<=0;
	end else begin
		dat1<=mdat1;
		dat2<=mdat2;
		dat3<=odat;
	end
end
wire [31:0] odat;
assign odat1 = odat;


fpmul mul0(
	.mdat1(dat1),
	.mdat2(dat2),
	.odat(odat),
	.clk(clk),
	.rst_n(rst_n)


);

/*
alt_fpmul mul1(
	.aclr(~rst_n),
	.clock(clk),
	.dataa(dat1),
	.datab(dat2),
	.result(odat),
	.zero());
*/
Endmodule