Перевод числа в строку с помощью SIMD + FPU

Вариант конвертирования double / real8 в строку посредством с SIMD использованием FPU в качестве вспомогательного средства. Применение FPU вызвано желанием получить 16 значащих цифр.

В соответствии с x64 software conventions будем считать что число подлежащие конвертированию расположено в XMM0.

Будем использовать x64 битный код при x32 битной адресации. Такой способ адресации позволяет использовать преимущества обоих диалектов.

Будем использовать недокументированное соглашение о передаче / возврате из функции множественных параметров. Соглашение абсолютно зеркально соглашению x64 software conventions за тем исключением что описывает правила размещения параметров при выходе из процедуры.

Для удобства чтения кода создадим два блока текстовых констант, в первом определим псевдонимы для аргументов ассемблерных команд, во втором псевдонимы размеров переменных в стеке которые позволят легко понимать какую именно переменную мы записываем / читаем:

ROUND_TOWARD_ZERO            equ 11b SIGNIFICANT_BIT_RESET        equ 3Fh

LCW          equ word LIExp2       equ dword LIExp10      equ dword LSExp10      equ dword LIUpPathNam  equ dword LILowPathNam equ dword LNamber      equ qword LMulExp2     equ qword LStX         equ tbyte LString      equ xmmword * 2

Создаем сегмент вспомогательных данных которые будем использовать при вычислении, особо стоит обратить внимание что данные для регистров SIMD выравнены по параграфу для возможности прямого обращения:

.data 	f10m4 real4 4 dup (1.0e-4) 	f10p4 real4 4 dup (1.0e+4) 	f10m2 real4 4 dup (1.0e-2) 	f10p2 real4 4 dup (1.0e+2) 	f10m1 real4 4 dup (1.0e-1) 	f10p1 real4 4 dup (1.0e+1)  	f0001 real4 0.0, 1.0e-2, 1.0e-1, 1.0 	f0002 real4 0.0, 0.0,    1.0e+1, 1.0e+1  	i30h  db 10h dup (30h)    	f10p8 real4 1.0e+8 	NoSD  real4 7.0  	CW0   dw 0F7Fh 	CW1   dw 037Fh 	DotM  dw 652Dh 	 	namber   real8 -1.234567890123456e+248

Для ускорения работы алгоритма отказываемся от полного сохранения среды FPU и ограничиваемся установкой расширенной точности, округлением к нулю и освобождением одного регистра:

	fstcw     word ptr[esp - LCW] 	fstp     tbyte ptr[esp - LCW - LStX] 	fldcw      CW0

Перегружаем Число в RAX и сбрасываем знаковый бит получая модуль Числа:

	movd      rax, xmm0 	btr       rax, SIGNIFICANT_BIT_RESET

Сохраняем в стек символ ‘-‘ предваряя его тремя нулевыми байтами, сохраняем в стек знаковый бит расширяя его до байта, а учитывая ранее загруженные три нулевых байта получаем в памяти двойное слов принимающие значение 1 для положительного числа и 0 для отрицательного и непосредственно сам модуль числа:

	mov      dword ptr[esp - LString - dword], 2D000000h ; '-' 00 00 00 	setnc     byte ptr[esp - LString - dword - byte]  	mov      qword ptr[esp - LCW - LStX - LNamber], rax

Извлекаем экспоненту числа и сохраняем ее в стек, особо стоит отметить что память в стеке всегда «горячая» и при «прямом» чтении / записи, то есть без «скачков», исключает ошибку кэш промаха, то есть чтение запись стабильно равны двум тактам, что конечно в два раза медленней чем работа с регистрами но при работе с FPU это не устранимая проблема:

	shr       rax, 34h 	sub       eax, 3FFh 	mov      dword ptr[esp - LCW - LStX - LNamber - LIExp2], eax

Загружаем в FPU логарифм 2 (двух) по основанию 10 (десять) и умножаем его на экспоненту Числа по основанию 2 (два), получая тем самым десятичный порядок Числа:

	fldlg2 	fimul    dword ptr[esp - LCW - LStX - LNamber - LIExp2]

Вычитаем из полученного порядка количество требуемых символов до запятой, и получаем десятичный порядок Множителя, который сохраняем в стек одновременно преобразовывая его в целое и округляя в сторону нуля:

	fsubr    NoSD 	fistp    dword ptr[esp - LCW - LStX - LNamber - LIExp10]

Загружаем в FPU логарифм 10 (десяти) по основанию 2 (два) и умножаем его на экспоненту Множителя по основанию 10 (десять), получая тем самым двоичный порядок Числа:

	fldl2t 	fimul    dword ptr[esp - LCW - LStX - LNamber - LIExp10]

Сохраняем в стек двоичный порядок Множителя одновременно преобразовывая его в целое и округляя в сторону нуля:

	fist     dword ptr[esp - LCW - LStX - LNamber - LIExp10 - LIExp2]

Вычитаем из двоичного порядка Множителя целую часть и находим двоичную мантиссу Множителя. Особо стоит обратить внимание что команда f2xm1 занимает до 60 тактов даже на Skylake и после нее разумно размещать код не требующий мгновенно результата команды f2xm1:

	fisub    dword ptr[esp - LCW - LStX - LNamber - LIExp10 - LIExp2] 	f2xm1

Загружаем в ХММ0 десятичный порядок Множителя одновременно преобразуя его в float:

	cvtsi2ss xmm0, dword ptr[esp - LCW - LStX - LNamber - LIExp10]

Загружаем в EAX двоичный порядок Множителя, находим его экспоненту и сохраняем ее в стек:

	mov       eax, dword ptr[esp - LCW - LStX - LNamber - LIExp10 - LIExp2] 	add        ax, 3FFh 	shl       rax, 34h 	mov            qword ptr[esp - LCW - LStX - LNamber - LMulExp2], rax

Сравниваем значение десятичной экспоненты числа с нулем и в случае равенства игнорируем участок кода по созданию строки экспоненты:

	xor       edx, edx 	subss    xmm0, NoSD 	pxor     xmm1, xmm1 	comiss   xmm1, xmm0 jz @f

Создаем вектор из четырех значений экспоненты и находим его модуль:

	shufps   xmm0, xmm0, 0 	subps    xmm1, xmm0 	maxps    xmm0, xmm1

Выделяем из трех значений экспоненты Числа три компоненты содержащие значение единиц десятков и сотен, а одно из значений экспоненты обнуляем за ненадобностью:

	mulps    xmm0, xmmword ptr f0001 	roundps  xmm0, xmm0, ROUND_TOWARD_ZERO 	pshufd   xmm1, xmm0, 10010000b 	mulps    xmm1, xmmword ptr f0002 	subps    xmm0, xmm1

Преобразуем компоненты экспоненты Числа в двойные слова и упаковываем их байты:

	cvtps2dq xmm0, xmm0 	pxor     xmm1, xmm1 	pcmpeqd  xmm1, xmm0 	packusdw xmm0, xmm0 	packuswb xmm0, xmm0

Загружаем в регистр AX строку 'е-' 'е+' в зависимости от знака экспоненты:

	mov       eax, 2B65h 	cmovc      ax, DotM

Создаем строку экспоненты Числа и помещаем ее в регистр RDX:

	movmskps  ecx, xmm1 	bsr       ecx, ecx 	lea       ecx,[ecx * 8 - 8] 	movd      edx, xmm0 	add       edx, 30303000h 	shrd      rdx, rdx, cl        	mov        dx,  ax

Умножаем мантиссу Множителя на экспоненту Множителя:

@@:	fmul     qword ptr[esp - LCW - LStX - LNamber - LMulExp2]

Добавляем к мантиссе Множителя экспоненту Множителя и получаем полное значение Множителя:

	fadd     qword ptr[esp - LCW - LStX - LNamber - LMulExp2]

Умножаем полное значение Множителя на модуль Числа и получаем число с восьми значащими числами до запятой и восьмью после:

	fmul     qword ptr[esp - LCW - LStX - LNamber]

Сохраняем в стек верхние восемь знаков целой части Числа одновременно преобразовывая его в целое и округляя в сторону нуля:

	fist     dword ptr[esp - LCW - LStX - LILowPathNam - LIUpPathNam]

Вычитаем из Числа целую часть и умножаем ее на 10е+8 перемещая нижние восемь чисел числа в целую часть Числа:

	fisub    dword ptr[esp - LCW - LStX - LILowPathNam - LIUpPathNam] 	fmul     f10p8

Устанавливаем округление к ближайшему числу и выгружаем нижние восемь чисел Числа одновременно преобразовывая его в целое и округляя в ближайшего числа:

	fldcw      CW1 	fistp    dword ptr[esp - LCW - LStX - LIUpPathNam]

Восстанавливаем среду FPU:

	fld      tbyte ptr[esp - LCW - LStX] 	fldcw     word ptr[esp - LCW]

Загружаем верхнюю и нижнюю часть Числа в регистр ХММ0 и конвертируем их в float:

	movq     xmm0, qword ptr[esp - LCW - LStX - LIUpPathNam - LILowPathNam] 	cvtdq2ps xmm0, xmm0

Разделяем верхнюю и нижнюю часть Числа на четыре части:

	movaps   xmm1, xmm0 	mulps    xmm0, xmmword ptr f10m4 	roundps  xmm0, xmm0, ROUND_TOWARD_ZERO 	movaps   xmm2, xmm0 	mulps    xmm2, xmmword ptr f10p4 	subps    xmm1, xmm2 	unpcklps xmm0, xmm1

Разделяем четыре части Числа на восемь частей:

	movaps   xmm1, xmm0 	mulps    xmm0, xmmword ptr f10m2 	roundps  xmm0, xmm0, ROUND_TOWARD_ZERO 	movaps   xmm2, xmm0 	mulps    xmm2, xmmword ptr f10p2 	subps    xmm1, xmm2

Разделяем нижние четыре части Числа на восемь чисел и преобразуем их в восемь целых слов:

	movaps   xmm2, xmm1 	mulps    xmm1, xmmword ptr f10m1 	roundps  xmm1, xmm1, ROUND_TOWARD_ZERO 	movaps   xmm3, xmm1 	mulps    xmm3, xmmword ptr f10p1 	subps    xmm2, xmm3	 	cvtps2dq xmm1, xmm1 	cvtps2dq xmm2, xmm2 	pslld    xmm2, 8 	paddb    xmm1, xmm2

Разделяем верхние четыре части Числа на восемь чисел и преобразуем их в восемь целых слов:

	movaps   xmm2, xmm0 	mulps    xmm0, xmmword ptr f10m1 	roundps  xmm0, xmm0, ROUND_TOWARD_ZERO 	movaps   xmm3, xmm0 	mulps    xmm3, xmmword ptr f10p1 	subps    xmm2, xmm3 	cvtps2dq xmm0, xmm0 	cvtps2dq xmm2, xmm2 	pslld    xmm2, 8 	paddb    xmm0, xmm2

Складываем верхнюю и нижнюю часть Числа:

	pslld    xmm1, 16 	paddb    xmm0, xmm1

Вычисляем длину строки в байтах:

	pxor     xmm3, xmm3 	pcmpeqb  xmm3, xmm0 	pmovmskb  eax, xmm3 	bts       eax, 10h 	bsr       eax, eax

Преобразуем числа в символы и сохраняем их в стек:

	paddb    xmm0, xmmWord ptr i30h 	movdqu   [esp - LString + byte], xmm0

Сохраняем строку экспоненты Числа в стек:

	mov      qword ptr[esp - LString + byte + eax], rdx

Вычисляем длину строки экспоненты Числа:

	movd     xmm0, rdx 	pxor     xmm1, xmm1 	pcmpeqb  xmm1, xmm0 	pmovmskb  edx, xmm1 	bsf       edx, edx 	lea       eax,[eax + edx + word + byte]

Вставляем символ ‘.’ между первым и вторым символом строки:

	mov     dl,[esp - LString + byte] 	mov     dh,'.' 	mov        [esp - LString], dx

Вычисляем полную длину строки:

	mov     ecx,   dword ptr[esp - LString - dword - byte] 	sub     eax, ecx

Сохраняем строку Числа в регистры ХММ1 и ХММ2:

	movdqu xmm1, xmmword ptr[esp - LString +   ecx - byte] 	movdqu xmm2, xmmword ptr[esp - LString +   ecx - byte + xmmword]

Дублируем значение длины Числа в регистр ECX:

	mov     ecx, eax

Зачем писать этот код если ранее ты уже разместил код про который заявил что он самый быстрый — потому что этот код еще быстрей.

Чем он лучше предыдущего — в этом коде векторизовано разложение Числа на числа.

Почему в нем одновременно используются FPU и SIMD — потому что в FPU есть режим расширенной точности позволяющий извлечь 16 значащих цифр.