倍精度2進浮動小数点数の10進変換 (2012/08/02)
Raspberry Pi は、普通にLinuxとしてできることはなんでもプログラミングできます。 Jun's homepage では Raspberry Pi のハードウェアに特化した部分から試してみたいと思います。候補は GPUを使った OpenGL ES2のプログラミングか、ARM11のコプロセッサであるVFPv2を使った浮動小数点演算をアセンブラで試すのが面白そうです。 VFPv2 はベクトル演算も可能なので色々実験するネタには困りません。
が、何をするにも浮動小数点演算の結果を表示して確認する必要があるので、最初のネタとして倍精度2進浮動小数点数を10進変換して標準出力に表示するサブルーチン(DoubleString)を作成してみました。 DoubleString ルーチンは d0 に設定された倍精度2進浮動小数点数を10進浮動小数点数の文字列に変換して、その文字列(AsciiZ)の先頭アドレスを r0 に返します。標準出力への表示のために、このページで作成したeabi 版の stdio.sを使用しています。
とにかく、Raspberry Pi での私の最初のプログラミングの成果として置いておきます。
@-------------------------------------------------------------------------
@ file : d2str.s 16digit version
@ 2012/08/02 Jun Mizutani (https://www.mztn.org/)
@-------------------------------------------------------------------------
.include "stdio.s"
.text
.global _start
_start:
adr r1, D001
1:
fldd d0, [r1] @ mem --> d0
add r1, r1, #8
bl DoubleString
bl OutAsciiZ
bl NewLine
fcmpezd d0 @ exit, if d0==0.0
fmstat
bne 1b @ next
bl Exit
@------------------------------------
D001: .double 1.2345678987654321
.double -1.2345678987654321
.double 9.87654E33
.double -9.87654E-33
.double 1.2345678987654321E299
.double -1.2345678987654321E-299
.double 12345678987654
.double 12345678987653
.double 12345678987652
.double -123456789876543
.double -123456789876542
.double -123456789876541
.double 923456789876543
.double 923456789876542
.double 923456789876541
.double 1.000
.double 0.1
.double 0
@-------------------------------------------------------------------------
@ float(d0) to asciiz string
@ r0 returns string address
@-------------------------------------------------------------------------
DoubleString:
stmfd sp!, {r1-r6, lr} @ push
vpush {d0-d5} @ fstmfdd sp!, {d1}
mov r4, #0 @ plus
fcmpezd d0 @ check 0
fmstat
beq 6f @ d0 = 0.0
fabsdlt d0, d0 @ if negative, d0 = abs(d0); r4=1
movlt r4, #1
ldr ip, numBuf
cmp r4, #0 @ if negative, print '-'s
mov r0, #' '
movne r0, #'-'
strb r0, [ip], #2 @ additional space for decimal point
@------------------------------
mov r6, #0 @ exponent
mov r3, #TblSize
mov r2, #512
fldd d3, DOne @ 1.00 --> d3
fldd d4, DTen @ 10.0 --> d4
fcmpd d0, d3
fmstat
bge 2f @ goto 2f, if d0 >= 1.0
1: @ negative exponent (f<1)
sub r3, r3, #8
mov r2, r2, ASR #1 @ start from w=256
cmp r2, #0
beq 4f @ range01_1
adr r5, NegTable
add r5, r5, r3
fldd d5, [r5]
fcmpd d0, d5 @ f > 10^-w
fmstat
sublt r6, r6, r2 @ exponent = exponent - w
fdivdlt d0, d0, d5 @ f = f / 10^-w = f * 10^w
b 1b
2: @ positive (f>=1)
sub r3, r3, #8
mov r2, r2, ASR #1 @ w
cmp r2, #0
beq 3f @ div 10
adr r5, PosTable
add r5, r5, r3
fldd d5, [r5]
fcmpd d0, d5 @ f <= 10^n
fmstat
addge r6, r6, r2 @ exponent = exponent + w
fdivdge d0, d0, d5 @ f = f / 10^w
b 2b
3: @ div10
fdivd d0, d0, d4 @ f = f / 10
add r6, r6, #1 @ w = w + 1
@------------------------------
4: @ mantissa should be 0.1 <= f < 1
bl VFPRoundTowardZero
fldd d2, Exp08
fmuld d0, d0, d2
ftosid s30, d0 @ d0 --> integer into s2
fmrs r4, s30 @ set upper 8 digit
fsitod d1, s30
fsubd d0, d0, d1 @ d0 = d0 - (upper 8)
fmuld d0, d0, d2 @ d0 = d0 * 10^8
bl VFPRoundToNearest
ftosid s30, d0 @ d0 --> integer into s30
fmrs r5, s30 @ set lower 8 digit
ldr r2, digit8
cmp r5, r2
addgt r4, r4, #1 @ add carry to upper 8
mov r0, r4
bl _put_uint @ upper 8 digits
mov r0, r5
bl _put_uint @ lower 8 digits
sub r0, r6, #1 @ r0 = r6-1
cmp r0, #0
movne r1, #'E'
strneb r1, [ip], #1
blne _put_int @ exponent
5:
ldr r0, numBuf
ldrb r1, [r0, +#2] @ shift left 1 digit for decimal point
strb r1, [r0, +#1]
mov r1, #'.'
strb r1, [r0, +#2]
vpop {d0-d5} @ fldmfdd sp!, {d1}
ldmfd sp!, {r1-r6, pc} @ pop & return
6: @ zero
ldr ip, numBuf
mov r0, #' '
strb r0, [ip], #2 @ additional space for decimal point
mov r0, #0
bl _put_uint @ upper 8 digits
bl _put_uint @ lower 8 digits
b 5b
@==============================================================
numBuf: .long num_buffer
digit8: .long 99999999
.align 3
DOne: .double 1.0
DTen: .double 10.0
Exp08: .double 1.0e8
NegTable:
.double 1.0e-1
.double 1.0e-2
.double 1.0e-4
.double 1.0e-8
.double 1.0e-16
.double 1.0e-32
.double 1.0e-64
.double 1.0e-128
.double 1.0e-256
PosTable:
.double 1.0e+1
.double 1.0e+2
.double 1.0e+4
.double 1.0e+8
.double 1.0e+16
.double 1.0e+32
.double 1.0e+64
.double 1.0e+128
.double 1.0e+256
TblSize = . - PosTable
@------------------------------------
@ Output 8 digit unsigned number to buffer
@ r0 : number
@ ip : buffer
_put_uint:
stmfd sp!, {r0-r3, lr} @ push
mov r3, #8 @
mov r2, #0 @ counter
1: mov r1, #10 @ r1 = 10
bl udiv @ division by 10
add r2, r2, #1 @ counter++
stmfd sp!, {r1} @ least digit (reminder)
cmp r2, r3
bne 1b @ done ?
2: ldmfd sp!, {r0} @ most digit
add r0, r0, #'0' @ ASCII
strb r0, [ip], #1 @ store a char into buffer
subs r2, r2, #1 @ counter--
bne 2b
strb r2, [ip] @ store 0 into buffer
ldmfd sp!, {r0-r3, pc} @ pop & return
@------------------------------------
@ Output signed number to buffer
@ r0 : number
@ ip : buffer
_put_int:
stmfd sp!, {r0-r2, lr} @ push
mov r2, #0 @ counter
cmp r0, #0
sublt r0, r2, r0
movlt r1, #'-'
strltb r1, [ip], #1
1: mov r1, #10 @ r1 = 10
bl udiv @ division by 10
add r2, r2, #1 @ counter++
stmfd sp!, {r1} @ least digit (reminder)
cmp r0, #0
bne 1b @ done ?
2: ldmfd sp!, {r0} @ most digit
add r0, r0, #'0' @ ASCII
strb r0, [ip], #1 @ store a char into buffer
subs r2, r2, #1 @ counter--
bne 2b
strb r2, [ip] @ store 0 into buffer
ldmfd sp!, {r0-r2, pc} @ pop & return
@-------------------------------------------------------------------------
@ Set VFP rounding mode
@-------------------------------------------------------------------------
VFPRoundTowardZero:
stmfd sp!, {r4, r5, lr}
mov r4, #0x00c00000 @ rounding mode(11)
fmrx r5, FPSCR @ fpscr --> r5
orr r5, r4 @ Round toward zero(RZ) mode
fmxr FPSCR, r5 @ r5 --> FPSCR
ldmfd sp!, {r4, r5, pc}
VFPRoundToNearest:
stmfd sp!, {r4, r5, lr}
mov r4, #0xff3fffff @ rounding mode(00)
fmrx r5, FPSCR @ fpscr --> r5
and r5, r4 @ Round to nearest(RN) mode
fmxr FPSCR, r5 @ r5 --> FPSCR
ldmfd sp!, {r4, r5, pc}
@==============================================================
.bss
.align 2
num_buffer: .skip 32
実行例
ソース中に埋め込んだ浮動小数点数が正しく変換されるか確認できます。十分な精度で変換できていると(私は)思います。
jun@raspberrypi ~/pi_asm/d2str $ make as -mcpu=arm1176jzf-s -mfpu=vfpv2 -o d2str.o d2str.s ld -o d2str d2str.o jun@raspberrypi ~/pi_asm/d2str $ ./d2str 1.234567898765432 -1.234567898765432 9.876540000000000E33 -9.876539999999999E-33 1.234567898765432E299 -1.234567898765432E-299 1.234567898765400E13 1.234567898765300E13 1.234567898765200E13 -1.234567898765430E14 -1.234567898765420E14 -1.234567898765410E14 9.234567898765430E14 9.234567898765421E14 9.234567898765408E14 1.000000000000000 1.000000000000000E-1 0.000000000000000
Makefile
アセンブルに使ったMakefileです。TABコードに注意。
# Makefile for vfp
# 2012/08/01 Jun Mizutani
OBJS = d2str.o
TARGET = d2str
ARCH=-mcpu=arm1176jzf-s -mfpu=vfpv2
AS=as
LD=ld
STRIP=strip
.s.o :
${AS} ${ARCH} -o $@ $<
all : ${TARGET}
${TARGET} : ${OBJS}
${LD} -o $@ ${OBJS}
strip : ${TARGET}
${STRIP} ${TARGET}
dist :
rm *.o
clean :
rm *.o