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

章國(guó)華

?

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

章國(guó)華

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

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

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

0 引言

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

假設(shè)：f1和f2是兩個(gè)浮點(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是否相等。如果要求更高的精度，也不能簡(jiǎn)單地把DBL_EPSILON定得更小就行了。

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

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

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

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

{

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

}

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

顯而易見，可以把除數(shù)選取為a和b當(dāng)中絕對(duì)數(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;

};

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

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

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

double f1 = 1.58;

double f2 = 1.57

根據(jù)以上的定義，事實(shí)上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ù)實(shí)際上是有限的，就是IEEE754枚舉的232個(gè)，事實(shí)上，絕大多數(shù)的浮點(diǎn)數(shù)數(shù)值是不能在計(jì)算機(jī)里精確表達(dá)的。

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

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

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

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

{

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

}

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

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

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

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

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

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

/* 輸出: TRUE，兩個(gè)浮點(diǎn)數(shù)64位近似相等；FALSE 兩個(gè)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;

}

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

#defineDBL_MIN

2.2250738585072014e-308

/*min positive value*/

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

3 結(jié)論

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

[1] 杜叔強(qiáng). 淺析C語言中的浮點(diǎn)數(shù)[J]. 蘭州工業(yè)高等專科學(xué)校學(xué)報(bào), 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ù)的存儲(chǔ)格式及其有效數(shù)字位數(shù)[J].計(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

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