BitVisorのcore/arith.sに含まれる各関数の実装をとても簡単に紹介します。
はじめに
/*
* Copyright (c) 2007, 2008 University of Tsukuba
* All rights reserved.
*
* Redistribution and use in source and binary forms, with or without
* modification, are permitted provided that the following conditions are met:
*
* 1. Redistributions of source code must retain the above copyright notice,
* this list of conditions and the following disclaimer.
* 2. Redistributions in binary form must reproduce the above copyright notice,
* this list of conditions and the following disclaimer in the documentation
* and/or other materials provided with the distribution.
* 3. Neither the name of the University of Tsukuba nor the names of its
* contributors may be used to endorse or promote products derived from
* this software without specific prior written permission.
*
* THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
* AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
* IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
* ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE
* LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
* CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
* SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
* INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
* CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
* ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
* POSSIBILITY OF SUCH DAMAGE.
*/
mpumul_64_64
64ビット同士の符号なし整数の掛け算をして128ビットの答えを出します。
void mpumul_64_64 (u64 m1, u64 m2, u64 ans[2]);
AMD64用
mpumul_64_64:
mov %rdi,%rax # m1 -> rax
mov %rdx,%rdi # ans -> rdi
mul %rsi # rax * m2 -> rdx:rax
mov %rax,0(%rdi) # rax -> ans[0]
mov %rdx,8(%rdi) # rdx -> ans[1]
ret
mul命令だけでできるので非常に簡単です。
IA-32用
mpumul_64_64:
mov 20(%esp),%ecx # ans
mov 4(%esp),%eax # m1l
mov 12(%esp),%edx # m2l
mul %edx # m1l * m2l
mov %eax,0(%ecx)
mov %edx,4(%ecx)
mov 8(%esp),%eax # m1h
mov 16(%esp),%edx # m2h
mul %edx # m1h * m2h
mov %eax,8(%ecx)
mov %edx,12(%ecx)
mov 4(%esp),%eax # m1l
mov 16(%esp),%edx # m2h
mul %edx # m1l * m2h
add %eax,4(%ecx)
adc %edx,8(%ecx)
adcl $0,12(%ecx)
mov 8(%esp),%eax # m1h
mov 12(%esp),%edx # m2l
mul %edx # m1h * m2l
add %eax,4(%ecx)
adc %edx,8(%ecx)
adcl $0,12(%ecx)
ret
mul命令を4回使っているのが見えます。やっていることは結構簡単です。64ビットのm1とm2を、32ビットずつのm1h, m1lおよびm2h, m2lに分けて、式を変形していきます。
m1*m2=(m1h*2^32+m1l)*(m2h*2^32+m2l)
=m1h*2^32*m2h*2^32+m1l*m2h*2^32+m1h*2^32*m2l+m1l*m2l
=m1h*m2h*2^64+(m1l*m2h+m1h*m2l)*2^32+m1l*m2l
となります。プログラム的には2のべき乗の掛け算がビットシフトになるので以下のような順番でやっています。
m1l*m2lを計算してans[0]に代入
m1h*m2hを計算してans[1]に代入
m1l*m2hを計算してans[0.5]に加算
m1h*m2lを計算してans[0.5]に加算
mpudiv_128_32
128ビットの符号なし整数を32ビットの符号なし整数で割って、商と余りを出します。
u32 mpudiv_128_32 (u64 d1[2], u32 d2, u64 quotient[2]);
AMD64用
mpudiv_128_32:
mov %rdx,%rcx # quotient -> rcx
xor %rdx,%rdx # 0 -> rdx
mov %esi,%esi # d2 (32bit) -> rsi
mov 8(%rdi),%rax # d1[1] -> rax
div %rsi # rdx:rax % rsi -> rdx, rdx:rax / rsi -> rax
mov %rax,8(%rcx) # rax -> quotient[1]
mov 0(%rdi),%rax # d1[0] -> rax
div %rsi # rdx:rax % rsi -> rdx, rdx:rax / rsi -> rax
mov %rax,0(%rcx) # rax -> quotient[0]
mov %rdx,%rax # return rdx
ret
2回div命令を使っています。これは、以下の筆算と同じような考え方によります。
16
-----
2 ) 33
2
----
13
12
----
1
この筆算では、33を2で割るのに、まず1桁の3を2で割り、1余り1の答えを出します。次に下の桁を入れて、13を2で割り、6余り1を出します。
10進数で1桁、0から10-1を取っているわけですが、これを2^64進数に当てはめると、1桁は0から2^64-1になります。まず上の64ビットを割ると上の答えが出て、次に余りを64ビット左にシフトして下の64ビットを加えて割ると、下の答えと余りが出るわけです。そういうプログラムになっているのがわかると思います。div命令というのは、%rdxを割られる数の上位ビットとしていながら、%rdxに余りを出すので、ビットシフトすら不要であり、最初からこういう使い方が想定されていたに違いありません。
IA-32用
mpudiv_128_32:
push %edi
push %esi
# 8=ret 12=d1[] 16=d2 20=quotient[]
mov 12(%esp),%edi # d1 -> edi
mov 16(%esp),%esi # d2 -> esi
mov 20(%esp),%ecx # quotient -> ecx
xor %edx,%edx # 0 -> edx
mov 12(%edi),%eax # ((u32*)d1)[3] -> eax
div %esi # edx:eax % esi -> edx, edx:eax / esi -> eax
mov %eax,12(%ecx) # eax -> ((u32*)quotient)[3]
mov 8(%edi),%eax # ((u32*)d1)[2] -> eax
div %esi # edx:eax % esi -> edx, edx:eax / esi -> eax
mov %eax,8(%ecx) # eax -> ((u32*)quotient)[2]
mov 4(%edi),%eax # ((u32*)d1)[1] -> eax
div %esi # edx:eax % esi -> edx, edx:eax / esi -> eax
mov %eax,4(%ecx) # eax -> ((u32*)quotient)[1]
mov 0(%edi),%eax # ((u32*)d1)[0] -> eax
div %esi # edx:eax % esi -> edx, edx:eax / esi -> eax
mov %eax,0(%ecx) # eax -> ((u32*)quotient)[0]
mov %edx,%eax # return edx
pop %esi
pop %edi
ret
なんか複雑そうに見えますが考え方はAMD64用と全く同じです。この関数はわざと割る数を32ビットにしてありますので、div命令を倍の4回使うだけで済みます。
ipchecksum
IPヘッダーのチェックサムの計算をします。1の補数による足し算です。
u16 ipchecksum (void *buf, u32 len);
AMD64用
ipchecksum:
mov %esi,%ecx # len (32bit) -> rcx
mov %rdi,%rsi # buf -> rsi
mov $-1,%rdi
xor %rdx,%rdx
cld
1:
shr %ecx # len bit0 test
jnc 1f
shl $8,%rdi
1:
test $6,%esi # rsi bit1 and bit2 test
je 1f
test %ecx,%ecx
je 1f
xor %eax,%eax
2:
lodsw
add %eax,%edx
sub $1,%ecx
je 1f
test $6,%esi
jne 2b
1:
shr %ecx # len bit1 test
jnc 1f
shl $16,%rdi
1:
shr %ecx # len bit2 test
jnc 1f
shl $32,%rdi
1:
not %rdi
test %ecx,%ecx
je 1f
2:
lodsq
add %rax,%rdx
adc $0,%rdx
sub $1,%ecx
jne 2b
1:
lodsq
and %rdi,%rax
add %rdx,%rax
adc $0,%rax
mov %eax,%edx
shr $32,%rax
add %edx,%eax
adc $0,%eax
mov %eax,%edx
shr $16,%eax
add %dx,%ax
adc $0,%ax
1:
xor $~0,%ax
je 1b
ret
なんだこれ... とにかく下のほうにあるlodsq, add, adcのところが基本です。なるべく高速化しようと思ってできるだけ8バイトずつ足すようにしています。あふれた桁がキャリーフラグに入り、そのフラグをadc命令で直接足すことができるので、アセンブリ言語で書くとシンプルです。その下にもlodsqがありますが、これは最後のはみ出した部分を処理するところです。読みすぎておいてビットマスクでデータを落として足し算をするみたいですね。そして32ビットに縮めた後16ビットに縮め、最後にビット反転しますが、0になったら0xffffにするように条件ジャンプが入っています。
IA-32用
ipchecksum:
push %edi
push %esi
mov 16(%esp),%ecx # len -> ecx
mov 12(%esp),%esi # buf -> esi
mov $-1,%edi
xor %edx,%edx
cld
1:
shr %ecx # len bit0 test
jnc 1f
shl $8,%edi
1:
test $2,%esi # esi bit1 test
je 1f
test %ecx,%ecx
je 1f
xor %eax,%eax
2:
lodsw
add %eax,%edx
sub $1,%ecx
je 1f
test $2,%esi
jne 2b
1:
shr %ecx # len bit1 test
jnc 1f
shl $16,%edi
1:
not %edi
test %ecx,%ecx
je 1f
2:
lodsl
add %eax,%edx
adc $0,%edx
jne 2b
1:
lodsl
and %edi,%eax
add %edx,%eax
adc $0,%eax
mov %eax,%edx
shr $16,%eax
add %dx,%ax
adc $0,%ax
1:
xor $~0,%ax
je 1b
pop %esi
pop %edi
ret
これも基本的にはAMD64用とやっていることは同じですが、32ビット単位 (lodsl) になるため短くなっています。
crc32
CRC-32の計算をします。
u32 crc32 (void *buf, u32 len);
AMD64用
crc32:
xchg %rsi,%rdi
xor %eax,%eax
test %edi,%edi
jne crc32_common
ret
ほとんどが共通処理になっています。長さ0バイトの時だけ何もせずにretします。
IA-32用
crc32:
push %edi
push %esi
xor %eax,%eax
mov 16(%esp),%edi # len -> edi
mov 12(%esp),%esi # buf -> esi
test %edi,%edi
je 1f
call crc32_common
1:
pop %esi
pop %edi
ret
こちらもスタックに積まれていた引数をレジスターに取り出すだけで、残りは共通処理です。
共通処理
# 32bit/64bit comon routine
# input: eax=0 esi=buf edi=len len>0
.align 16
crc32_common:
cld
not %eax
3:
mov %eax,%edx
lodsb
mov $8,%cl
xor %al,%dl
movzbl %dl,%eax
shr %cl,%edx
2:
shr %eax
jnc 1f
xor $0xEDB88320,%eax
1:
sub $1,%cl
jne 2b
xor %edx,%eax
sub $1,%edi
jne 3b
not %eax
ret
あぁ、コメントのスペルが間違っている! ('A`)
1バイトずつ読み出して、ビットシフトを用いて1ビットずつ調べながらxorを行っています。事前にテーブル作成等を行わないとても遅いルーチンになっています。