在Ｃ語言中雙精度浮點(diǎn)數(shù)線性化相等比較的研究

章國華

?

在Ｃ語言中雙精度浮點(diǎn)數(shù)線性化相等比較的研究

章國華

（武漢船舶職業(yè)技術(shù)學(xué)院，武漢430050 ）

以LCC軟件為設(shè)計基礎(chǔ)，提出了雙精度浮點(diǎn)數(shù)比較大小算法的設(shè)計方法，通過程序驗證了算法的精確度，分別分析了絕對誤差和相對誤差方法設(shè)計算法，最后從雙精度浮點(diǎn)數(shù)在計算機(jī)中的表示方法角度提出了浮點(diǎn)數(shù)線性化后相等比較改進(jìn)的方法。

比較 算法 雙精度 有效數(shù)字

0 引言

Ｃ程序設(shè)計的教科書中一般都不涉及浮點(diǎn)數(shù)的運(yùn)算，浮點(diǎn)數(shù)的運(yùn)算涉及到計算機(jī)的復(fù)雜的硬件結(jié)構(gòu)和理解起來有一定難度的表示方法[1]。然而浮點(diǎn)數(shù)的運(yùn)算在計算機(jī)軟件系統(tǒng)中可以說是無處不在的。大多數(shù)編程語言都將浮點(diǎn)數(shù)據(jù)類型作為基本類型，從計算機(jī)的CPU硬件到軟件編譯器、操作系統(tǒng)都涉及到浮點(diǎn)數(shù)的運(yùn)算。當(dāng)然，在用各種編程語言實現(xiàn)各種應(yīng)用程序設(shè)計的時候，即數(shù)據(jù)的處理，大都會碰到浮躁數(shù)的運(yùn)算[2]。例如，在用C語言寫應(yīng)用軟件的時候，各類整數(shù)的相等比較是直接使用“==”來判斷，正是由于浮點(diǎn)數(shù)在計算機(jī)中表示與整數(shù)在計算機(jī)中的表示具有完全不一樣的特殊性，浮點(diǎn)數(shù)的比較就完全不能簡單的用“==”來實現(xiàn)。通常會利用差值的絕對值的精度來判斷。

假設(shè)：f1和f2是兩個浮點(diǎn)數(shù)，在Ｃ語言中規(guī)定有雙精度精度DBL_EPSILON。

#defineDBL_EPSILON

2.2204460492503131e-016

/*smallest such that

1.0+DBL_EPSILON!=1.0*/

即可以用 fabs(f1-f2)<= DBL_EPSILON 來判斷f1和f2是否相等。如果要求更高的精度，也不能簡單地把DBL_EPSILON定得更小就行了。

因為，DBL_EPSILON是C語言編譯器給出的一個不變的數(shù)據(jù)，也就是誤差分析當(dāng)中的絕對誤差，使用一個完全固定的數(shù)值，對于雙精度類型的數(shù)據(jù)可以表達(dá)的整個數(shù)域來說是有缺陷的。例如對于f1和f2大小分別是１e-100 、２e-100附近這樣的數(shù)據(jù)的時候，它是不合適的，因為它們之差1e-100已經(jīng)遠(yuǎn)小于DBL_EPSILON，兩個不等的浮點(diǎn)數(shù)就可以判斷是相等的，出現(xiàn)了相等比較的錯誤。即使DBL_EPSILON已經(jīng)是最小的，而且它是由數(shù)據(jù)在計算機(jī)中表示時的位數(shù)決定的。例如，對于f1和f2大小是10和10.0000000000000000001這樣的數(shù)據(jù)的時候，它也不合適，因為10和10.0000000000000000001兩個數(shù)的絕對誤差為1e-19,小于DBL_EPSILON，判斷的結(jié)果當(dāng)然是相等，也就出錯了，因為兩個數(shù)據(jù)的有效數(shù)字位數(shù)超過了16位。適合浮點(diǎn)數(shù)相等比較的情況只是f1或者f2在數(shù)值1附近的時候，函數(shù)的注釋正好說明了這個問題，大范圍的浮點(diǎn)數(shù)的比較需要解決。

1 雙精度浮點(diǎn)數(shù)比較相等方法的改進(jìn)

以上用絕對誤差的方式比較浮點(diǎn)數(shù)的相等不合適，就用相對誤差方式嘗試解決浮點(diǎn)數(shù)相等比較的問題。相對誤差的算法如下：

bool IsEqual(double a, double b, double relError )

{

return （( fabs ( (a-b)/a ) < relError )？TRUE ：FALSE;）

}

這個浮點(diǎn)數(shù)比較相等的函數(shù)也是有問題的。用固定的函數(shù)的第一個形參做相對比較，在應(yīng)用過程中，調(diào)用IsEqual(a, b, relError ) 和 IsEqual(b, a, relError ) 的時候，由于參數(shù)a和b的數(shù)量級不一樣，得到的結(jié)果是不同的。顯然，當(dāng)實際調(diào)用的第一個參數(shù)是0的話，就有出現(xiàn)了除0溢出的錯誤，這樣的相對誤差比較顯然也是不合適的。

顯而易見，可以把除數(shù)選取為a和b當(dāng)中絕對數(shù)值較大的就可以避免以上情況[3]：

bool IsEqual(double a, double b, relError )

{

if (fabs(a)

return ( fabs((a-b)/a) > relError ) ？TRUE ： FALSE;

else

return (fabs( (a-b)/b) > relError ) ？TRUE ： FALSE;

};

此算法就克服了絕對誤差算法對于太小的浮點(diǎn)數(shù)不能比較的缺陷。

2 優(yōu)化浮點(diǎn)數(shù)比較大小的實現(xiàn)方法

相對誤差算法能處理很小的浮點(diǎn)數(shù)相等的比較，但浮點(diǎn)數(shù)的運(yùn)算速度慢，能改用整數(shù)運(yùn)算就快了。為了從根本上解決問題，需要仔細(xì)研究浮點(diǎn)數(shù)在計算機(jī)中的表示。根據(jù)IEEE754浮點(diǎn)數(shù)的內(nèi)存結(jié)構(gòu)，符號1位在最高位，指數(shù)11位在次高位，尾數(shù)52位在低位（不包括隱含的１位）。浮點(diǎn)數(shù)的大小對應(yīng)的是其指數(shù)和尾數(shù)的大小，人工比較時，先比較符號，再比較指數(shù)，最后比較尾數(shù)，就可知兩個浮點(diǎn)數(shù)是否相等。用C語言編程來比較浮點(diǎn)數(shù)的相等也應(yīng)該是同樣的道理。即根據(jù)浮點(diǎn)數(shù)的內(nèi)存結(jié)構(gòu)，按照整數(shù)來理解進(jìn)行相等的比較，情況也是相同的。而且用這種方法對浮點(diǎn)數(shù)進(jìn)行比較的話，因為作為整數(shù)運(yùn)算效率比浮點(diǎn)數(shù)的高得多。比如

double f1 = 1.58;

double f2 = 1.57

根據(jù)以上的定義，事實上f1>f2 是成立，能否用一種方法，證明(long long int)f1 > (long long int)f2 也是成立的。根據(jù)IEEE754的浮點(diǎn)結(jié)構(gòu)特點(diǎn)，不是所有的浮點(diǎn)數(shù)都可以精確的表達(dá)，可以精確表達(dá)的浮點(diǎn)數(shù)實際上是有限的，就是IEEE754枚舉的232個，事實上，絕大多數(shù)的浮點(diǎn)數(shù)數(shù)值是不能在計算機(jī)里精確表達(dá)的。

目前雖然有80位的浮點(diǎn)數(shù)表示，這里只分析64位的情況（IEEE754）。尾數(shù)是53位的（暗含了第一個位數(shù)是1），對于可以精確表達(dá)的浮點(diǎn)數(shù)來說，對于IEEE754 單精度和雙精度浮點(diǎn)數(shù)，能夠精確表示的整數(shù)的范圍為[4]：

如果把這53位當(dāng)作整數(shù)來理解，對于隱含的１位，除浮點(diǎn)數(shù)零外，其它數(shù)都是１,不考慮也不影響比較的結(jié)果。對于任一浮點(diǎn)數(shù)，把它當(dāng)作64位的整數(shù)，根據(jù)二進(jìn)制計數(shù)原理，而且是線性分布的。這樣，將兩個浮點(diǎn)數(shù)看成是整數(shù)，不需要轉(zhuǎn)換成整數(shù)，將其對應(yīng)的整數(shù)做差值運(yùn)算，得到的整數(shù)表明的是兩個浮點(diǎn)數(shù)之間有多少個實際可以精確表達(dá)的浮點(diǎn)數(shù)的個數(shù)（對應(yīng)的指數(shù)相同的情況下是不言而喻的；指數(shù)不同的時候，也是同樣有效的，因為現(xiàn)在研究的是相等比較，不是大小比較，指數(shù)不同則涉及的是精確表達(dá)的浮點(diǎn)數(shù)的個數(shù)更多了，有利于相等的比較，而不會有相反的作用）。因此，對于兩個正的浮點(diǎn)數(shù)，他們的大小比較就可以用 (long long int)f1 - (long long int)f2 來進(jìn)行比較了，即將兩個浮點(diǎn)數(shù)當(dāng)成整型數(shù)來相減，其差值的結(jié)果實際上就相當(dāng)于相對誤差了，這個相對誤差，不等同于普通意義上的相對誤差，它所表達(dá)的是，兩個浮點(diǎn)數(shù)之間可能還有多少個可以精確表達(dá)的浮點(diǎn)數(shù)。這樣通過指定這個閾值來控制兩個浮點(diǎn)數(shù)的比較就更有效了[5]。(long int)f1 - (long int)f2在這里只是算法的示意，Ｃ編譯器會報錯。具體實現(xiàn)見最后面的程序。

對于兩個正的浮點(diǎn)數(shù)比較相等就有：

bool IsEqual(double f1, double f2, int precisionFloatNum)

{

if ( abs ( (long long int)f1 - (long longint)f2 ) < precisionFloatNum ) return TRUE;

}

這里要用整數(shù)的絕對值比較函數(shù)而不能用浮點(diǎn)數(shù)的絕對值比較函數(shù)，因為要將此時的浮點(diǎn)數(shù)看成是整數(shù)，否則就無法實現(xiàn)以上分析的目的。但是負(fù)整數(shù)和正整數(shù)之間現(xiàn)在還不能進(jìn)行直接的比較，因為根據(jù)IEEE754的內(nèi)存結(jié)構(gòu)，正整數(shù)和負(fù)整數(shù)是不同的，對應(yīng)的在計算機(jī)中表示的整數(shù)在整個區(qū)間是不連續(xù)的，正的最小的整數(shù)就是0，對應(yīng)的計算機(jī)表示的整數(shù)是0x0000000000000000，負(fù)的最大的整數(shù)就是－0，與之對應(yīng)的計算機(jī)表示的整數(shù)則是0x 8000000000000000，在IEEE754的表達(dá)當(dāng)中是有兩個0的，一個是 +0 一個是-0。而且按照 f1 == f2 的判斷 +0和-0是相等的，根據(jù)浮點(diǎn)數(shù)對應(yīng)的在計算機(jī)中表示的整數(shù)在整個區(qū)間的不連續(xù)性，具有以下特點(diǎn)，+0 和正的計算機(jī)表示的浮點(diǎn)數(shù)可以當(dāng)成整數(shù)的方式直接進(jìn)行比較，-0 和負(fù)的計算機(jī)表示的浮點(diǎn)數(shù)可以當(dāng)成為整數(shù)的方式直接進(jìn)行比較，如果采取措施將把他們統(tǒng)一起來，計算機(jī)表示的浮點(diǎn)數(shù)可以當(dāng)成整數(shù)直接進(jìn)行整數(shù)比較方式算法就完備了。仔細(xì)分析負(fù)整數(shù)的結(jié)構(gòu)，把負(fù)整數(shù)對應(yīng)的整數(shù)減去 -0 ，這樣，浮點(diǎn)數(shù)對應(yīng)的在計算機(jī)中表示的整數(shù)在整個區(qū)間就連續(xù)了。因此，所有的計算機(jī)表示用補(bǔ)碼的負(fù)整數(shù)經(jīng)過這次減法后，對應(yīng)的整數(shù)也都是原碼的負(fù)整數(shù)了，這樣整個整數(shù)比較就變得連續(xù)了就相當(dāng)于在整個浮點(diǎn)數(shù)范圍內(nèi)都是有效的了，只有被比較的浮點(diǎn)數(shù)連續(xù)或線性分布，將浮點(diǎn)數(shù)看成整數(shù)的比較算法才有意義。算法的關(guān)鍵是負(fù)數(shù)只需將最高位的１去掉，方法見程序的實現(xiàn)。注意，這里只是比較相等，不是比較大小，將被比較的兩個64位浮點(diǎn)數(shù)之間允許有多少個可以精確表達(dá)的浮點(diǎn)數(shù)作為比較相等的關(guān)鍵，從而實現(xiàn)快速和高精度的算法。

最后的浮點(diǎn)數(shù)比較算法就是：

/* 函數(shù): bool IsEqual(double f1, double f2, int precisionFloatNum) */

/* 功能：兩個64位浮點(diǎn)數(shù)是否近似相等*/

/* 輸入：兩個64位浮點(diǎn)數(shù)f1, f2 */

/* precisionFloatNum 被比較的兩個64位浮點(diǎn)數(shù)之間允許有多少個可以精確表達(dá)的浮點(diǎn)數(shù) */

/* 輸出: TRUE，兩個浮點(diǎn)數(shù)64位近似相等；FALSE 兩個64位浮點(diǎn)數(shù)不等 */

bool IsEqual(double f1, double f2, long int precisionFloatNum)

{

void *p1=&f1;

void *p2=&f2;

long long int bits1=*((long long int*)p1);

long long int bits2=*((long long int*)p2);

if((bits1>0&&bits2>0)||(bits1<0&&bits2<0))

{ bits1=(bits1>0)?bits1:(bits1-0x8000000000000000);

bits2=(bits2>0)?bits2:(bits2-0x8000000000000000);

return ((abs(bits1-bits2))< precisionFloatNum) ? TRUE : FALSE;

}

else

return FALSE;

}

通過編程實現(xiàn)浮點(diǎn)數(shù)相等的比較，相對誤差比較的算法和將浮點(diǎn)數(shù)當(dāng)在整數(shù)即線性化的比較的算法都能實現(xiàn)16位有效數(shù)字的精度，可實現(xiàn)最小到DBL_MIN數(shù)的比較。DBL_MIN的定義是：

#defineDBL_MIN

2.2250738585072014e-308

/*min positive value*/

在浮點(diǎn)數(shù)的絕對值大于DBL_MIN的基礎(chǔ)上，所有浮點(diǎn)數(shù)相等的比較現(xiàn)在只取決于有效數(shù)的位數(shù)，由此完整實現(xiàn)了浮點(diǎn)數(shù)的相等比較。

3 結(jié)論

本文從浮點(diǎn)數(shù)在計算機(jī)內(nèi)部表示的角度出發(fā)，介紹了雙精度浮點(diǎn)數(shù)比較大小的優(yōu)化算法，通過程序調(diào)試驗證了算法的有效性，并提出了改進(jìn)方法。對于兩個數(shù)據(jù)的有效數(shù)字位數(shù)超過了16位(253)的浮點(diǎn)數(shù)的相等比較，只有提高數(shù)據(jù)表示的位數(shù)，才能實現(xiàn)，例如，長雙精度浮點(diǎn)數(shù)有80位的數(shù)據(jù)寬度。從運(yùn)算速度上看，將浮點(diǎn)數(shù)當(dāng)成整數(shù)的浮點(diǎn)數(shù)線性化的算法在運(yùn)算速度和精度上優(yōu)于絕對和相對精度的算法，因為整數(shù)的減法和取絕對值運(yùn)算速度比浮點(diǎn)數(shù)的運(yùn)算要快。這里用LCC編譯器實現(xiàn)的雙精度的浮點(diǎn)數(shù)比較相等，其它編譯環(huán)境的算法則需要一定的修改。至于單精度的浮點(diǎn)數(shù)和長雙精度的浮點(diǎn)數(shù)，修改程序的相關(guān)部分即可實現(xiàn)。

[1] 杜叔強(qiáng). 淺析C語言中的浮點(diǎn)數(shù)[J]. 蘭州工業(yè)高等專科學(xué)校學(xué)報, 2010,17(5):27.

[2] David Goldberg. What every computer scientist should know about doubling-point Arithmetic [M]. Computing Surveys, 1991.

[3] 張宗杰,張明亮. C 語言中浮點(diǎn)數(shù)的存儲格式及其有效數(shù)字位數(shù)[J].計算機(jī)與數(shù)字工程, 2006,34(1):85.

[4] http://blog.csdn.net/seizef/article/details/5571783.

[5] http://cmdblock.blog.51cto.com/415170/600378.

[6] 陳煒峰. 電磁脈沖模擬器及其應(yīng)用研究[J]. 南京: 東南大學(xué)博士學(xué)位論文, 2007.

Approach to Compare to Equation of Double-precision Floating Point Numbers of Linearization in C Language

Zhang Guohua

(Wuhan Institute of Shipbuilding Technology, Wuhan430050, China)

Based on LCC software for design, this paper puts forward the design method of the double precision doubling- point comparisons algorithm, the accuracy of the algorithm is verified in applications, and the algorithm is designed with the absolute error and relative error method respectively. Finally aimed at representation for double-precision doubling-point number, the paper presents the arithmetic for floating point linearization.

compare; algorithms; double precision; significance digit

TP 312.1.4

A

1003-4862（2017）01-0040-03

2016-11-15

章國華(1964-),男，副教授。研究方向：機(jī)電一體化技術(shù)教學(xué)與研究。Email:993468391@qq.com