Pythonで小数を2進数に変換する
Pythonで小数点の計算をすると、おかしなことが起こります。極めて小さな誤差なので多くの場合は大した問題ではありませんが、やはり気になります。そこで、小数点以下の数値の取り扱いについて明らかにしていきます。
Pythonで小数点
0.1を単純に60桁まで表示させると、小数点第18位から555111・・・とごくわずかですが誤差が生じてしまいます。これは、Python、というよりコンピュータが2進数で数値を計算しており、そこでは浮動小数点方式という計算方法によっているためです。
- #1 0.1の小数表示
- print(f'{0.1:.60f}')
- #0.100000000000000005551115123125782702118158340454101562500000
そこで、なぜ0.1がこのようになってしまうのか、2進数に変換して確認してみます。
倍精度浮動小数点方式
ここでは簡単のため、整数の部分が0である小数、例えば0.625を2進数に変換します。
小数を2進数に変換する
整数の場合には、$2^1$、$2^2\cdots$の桁に”1”が立つかを計算していきます。これに対し、小数を2進数に変換する場合、$2^{-1}=0.5$、$2^{-2}=0.25$、$2^{-3}=0.125$のそれぞれの桁に”1”が立つかどうかで計算していきます。
はじめに、0.625を$2^{-1}=0.5$で割り商(いくつ含まれているか)を求め、その値を小数点第1位とします。2進数の商は1か0なので、0.625が$2^{-1}=0.5$と比べて大きいか等しければ”1”になるとも考えられます。小数点第1位で”1”が立つので、小数点第2位を求めるため$2^{-1}=0.5$との差、つまり0.125を計算します。
次に、0.125と$2^{-2}=0.25$との大小(商が1か?)を判断し、小数点第2位を求めます。ここでは、0.125は$2^{-2}=0.25$より小さいので”0”になります。さらに0.125を$2^{-3}=0.125$と比較すると等しくなるので(商が1)、小数点第3位は”1”になります。ここで余りは0なので、割り切れたと判断し計算は終了です。結果として、0.625は$2^{-1}=0.5$の桁と$2^{-3}=0.5$の桁が1、$2^{-2}=0.25$の桁が0となるので0.101となります。
上記の考え方を踏まえ、次の流れで小数を2進数に変換するとプログラミングしやすくなります。
Pythonで小数を2進数に変換する
前節での考え方をもとに、10進数の小数を2進数に変換する関数を作成します。なお、ここでは有効桁数(はじめに”1”が立ってからの桁数)を55桁として計算します。55桁というのは倍精度浮動小数点の仮数部が53桁なので、念のためこれより少し多めという意味です。
- #2 10進数の小数点以下の数値を2進数に変換
- from decimal import Decimal
- def bin_f(f):
- dec=Decimal(str(f))
- Bin=''
- nth=0
- first=0
- while dec:
- if dec >= 1:
- Bin += '1'
- dec = dec -1
- if first==0:
- first=nth
- else:
- if nth!=0:
- Bin += '0'
- else:
- Bin +='0.'
- dec*=2
- nth+=1
- if nth-first==55:
- break
- return Bin
- print(bin_f(0.625))
- print(bin_f(0.1))
0.101 0.0001100110011001100110011001100110011001100110011001100110
4. 引数として受け取った数値は、誤差が生じないようDecimal関数を使いdecとします。
5. 2進数に変換後の値を計算するため、変数Binを文字列として定義します。
6. 小数点第何位を計算しているのかを把握するため、変数nthを整数型として定義します。
7. はじめに”1”になるのが何桁目かを記録しておくため、変数firstを整数型として定義しておきます。
8. while文で、decが0になり割り切れるまで計算を繰り返します。なお、while文の中で順次decを2倍し、値が1以上であればBinの末尾に”1”を追加するとともに1を差し引き、1未満であれば”0”を追加する処理を繰り返します。
21. 割り切れないときは無限ループになることを防ぐため、nth-firstで初めに”1”が立ってからの桁数を計算し、この桁数が55になったらプログラムを終了するようにします。
浮動小数点型の変数の定義
前節で計算した2進数への変換のプログラムを変更し、正規化、つまり仮数部(Bin)と指数部(first)を分けて計算するようにします。このため、前のプログラムの仮数部の計算を、はじめて”1”が立ってから足しこんでいくように変更します。
- #3 10進数の小数点以下の数値を2進数に変換
- from decimal import Decimal
- def bin_f2(f):
- dec=Decimal(str(f))
- Bin=''
- nth=0
- first=0
- resudial=False
- while dec:
- if dec >= 1:
- Bin += '1'
- dec = dec -1
- if first==0:
- first=nth
- else:
- if first !=0:
- Bin += '0'
- dec*=2
- nth+=1
- if nth-first==55:
- break
- return Bin,first
- print(bin_f2(0.625))
- print(bin_f2(0.1))
('101', 1) ('1100110011001100110011001100110011001100110011001100110', 4)
7. 一番初めに”1”になる桁数を記録するため、変数firstを整数0として定義します。計算の途中である桁が”1”になったとき、firstが0であれば初めて”1”が立ったと判断し、そのときの桁数(nth)を代入します。firstは最終的にはバイアスの桁数となります。
16. 正規化の考え方から、はじめて”1”になった桁が仮数部の先頭となるため、firstが0のうちは”0”がBinに追加されないようにします。
22. 上記の変更により、仮数部と指数部の組合せをタプルで返ることができます。
さらに2進数を計算するプログラムを改良する
小数を2進数に変換する際に、分数1/3と同じように永久に割り切れないことが多くあります。このような場合に対応できるよう、プログラムを改良します。
改良されたプログラム
前節のプログラムに対し、53桁まで計算しても割り切れない場合には丸め処理をするようにします。また、0.625のように早々に割り切れた場合には、53桁目まで”0”を付け加えるように改良します。このことにより、計算結果とPythonで内部的に把握している値とを比較することができます。さらに仮数部の桁数を53桁だけではなく、自由に設定できるようにします。
- #4 10進数の小数点以下の数値を2進数に変換 ( 丸め 0をつける)
- from decimal import Decimal
- def bin_f3(f,degits=53):
- dec=Decimal(str(f))
- Bin=''
- nth=0
- first=0
- resudial=False
- while dec:
- if dec >= 1:
- Bin += '1'
- dec = dec -1
- if first==0:
- first=nth
- else:
- if first !=0:
- Bin += '0'
- dec*=2
- nth+=1
- if nth-first==degits+1:
- if dec == 0:
- resudial=False
- else:
- resudial=True
- Bin=round_bin(Bin,resudial)
- break
- return Bin+'0'*(degits-len(Bin)),first
- print(bin_f3(0.625))
- print(bin_f3(0.1))
('10100000000000000000000000000000000000000000000000000', 1) ('11001100110011001100110011001100110011001100110011010', 4)
前節のプログラムからの変更点は次の通りです。
3. 引数として、数値の値(f)と仮数の桁数(digits)を受け取ります。
20. 求めたい桁数(digits)より1桁多く計算しても割り切れないときは、round_rn関数を使い丸め処理をします。round_rn関数は次の節に記載します。
28. 計算結果として仮数部と指数部の桁数を返します。このとき、仮数部(Bin)が53桁に達する前に割り切れたときは53桁まで"0”を補います。
丸めを計算する関数
前節のプログラムで使用する丸め処理をする関数です。小数を2進数に変換するときに、倍精度浮動小数点方式では53bitを仮数部として使いますが、ここまで計算しても割り切れない場合には、最後のbitについて丸め処理をします。
丸め処理において、仮数部(有効な桁)の最後の桁をulp(Unit in the Last Place)といいます。また、その次のbitをgb(Guard bit)といい、さらにそれ以下のbitのうち1つでも”1”が立てば”1”になるようなbitをrb(Round bit)といいます。実際に計算するときはgbまで計算し、ここで割り切れればrb=”0”、割り切れなければrb=”1”とします。
丸め処理ではgbとrbを組み合わせ、その数値が切り上げた値と切り捨てた値のどちらに近いかで判断します。まず、gb=”1”、rb=”1” のときには切り上げた値に近いと考えます。ところが、rb=”0”、つまりgbの桁で割り切れたときには引き分け(どちらとも同じ近さ)となります。ここで常に切り上げるなどと決めてしまうと、全体として数値が大きい方に偏ってしまいます。そこで、ulpが”1”のときには切り上げ、ulpが”0”のときには切り捨てとして、uplが”0”になるように丸めます。これは、10進数で最後の桁が偶数になるように丸めるのと同じで、2進数の場合も、最後の桁がなるべく0になるようにしていると考えられます。一方、gb=”0”の場合は切り捨てた値に近づくので、無条件に切り捨てとなります。これらを表にすると、次の通りです。
ulp | gb | rb | 結果 |
---|---|---|---|
0 | 0 | 0 | 切り捨て |
0 | 0 | 1 | 切り捨て |
0 | 1 | 0 | 切り捨て |
0 | 1 | 1 | 切り上げ |
1 | 0 | 0 | 切り捨て |
1 | 0 | 1 | 切り捨て |
1 | 1 | 0 | 切り上げ |
1 | 1 | 1 | 切り上げ |
- def round_bin(Bin,rb): #1桁で判断 決定版
- ulp=Bin[-2] #unit in the last place
- gb =Bin[-1] #Guard bit
- sw_up= False
- if gb=='1':
- if rb==True: #Round bit
- sw_up=True
- else:
- if ulp=='1':
- sw_up=True
- result=Bin[:-1]
- if sw_up==True:
- return bin(eval('0b'+result) +eval('0b'+(len(Bin)-2)*'0'+'1'))[2:]
- else:
- return result
- print(round_bin('110',True))
- print(round_bin('111',True))
- print(round_bin('111',False))
- print(round_bin('101',False))
- print(round_bin('110',False))
- print(round_bin('100',False))
101 101 110 110 101
2. 引数Binで仮数部+gb(求めたいbitより1桁だけ多く)までの2進数を文字列として、またgbまで計算して割り切れなかった場合には、rb=Trueとして受け取ります。
5. sw_upはいったんFlaseで定義し、切り上げると判断したときにはTrueを代入するようにします。このため、gb=”0”の場合のように、切り上げの条件にあてはまらなければ自動的に切り捨てられるようになります。
6. gb=”1”のとき、rb=Trueであれば1つ繰り上がるので、sw_up=Trueとします。
9. gb=”1”でrb=Falseのときは、ulp=”0”になるように丸めるため、ulp=”1”の場合に切り上げるようにします。
11. resultにはBinからgbを取り除いた値が代入されます。したがって切り捨てられるときには、この値が戻り値になります。
12. 切り上げられる場合は、resultを2進数に変換した数値と、仮数部で最小の値(000・・・01)を2進数に変換した数値の和を計算します。この値が切り上げられるときの戻り値になります。
2進数に変換する関数を検証する
小数を2進数に変換するbin_3関数が正しいのか検証してみます。まず、Pythonのhexメソッドを通して2進数に変換するhex2bin関数を作成します。
こちらは少なくとも16進数に変換するまでは正しく、2進数に変換する部分はそんなに複雑なことはやっていないので正しく変換できているものと思われます。
- #5 16進数から2進数に変換
- def hex2bin(h):
- h=h.hex()
- Bin='1'
- for i in range(4,17):
- Bin+=(format(int(h[i],16),'#06b')[2:6])
- return Bin,int(h[19:])
- print(hex2bin(0.625) )
- print(hex2bin(0.1) )
('10100000000000000000000000000000000000000000000000000', 1) ('11001100110011001100110011001100110011001100110011010', 4)
ここまでは合ってそうです。
iに100万まで代入し1を割ることにより、100万分の1ごとにbin_f3関数とhex2bin関数を計算し、値が一致しないときには出力することにします。
- #6 確認
- for i in range(1,1000000):
- x=i/1000000
- if bin_f3(x)!=hex2bin(x):
- print(x,bin_f3(x),hex2bin(x))