From Wiki.conus.info

FPU (Floating-point unit) - девайс в процессоре работающий с числами с плавающей запятой. Раньше он назывался сопроцессором. Он немного похож на программируемый калькулятор и стоит немного в стороне от основного процессора. Интересен факт, что в свое время (до 80486) сопроцессор был отдельным чипом и на материнской плате он стоял не всегда. Раньше, его можно было докупить отдельно и поставить.

Перед изучением FPU полезно ознакомиться с тем как работают стековые машины:

http://en.wikipedia.org/wiki/Stack_machine

Или ознакомиться с основами языка Forth:

http://en.wikipedia.org/wiki/Forth_(programming_language)

Начиная с 80486, FPU уже всегда входит в состав процессора.

FPU имеет стек из восьми 80-битных регистров, каждый может содержать число в формате IEEE 754.

В Си типы имеются два типа для работы с числами с плавающей запятой, это float (32 бита) и double (64 бита). GCC поддерживает тип long double (80 бит), но MSVC - нет. Не смотря на то что float занимает столько же места сколько int на 32-битной архитектуре, представление чисел совершенно другое. Число состоит из знака, мантиссы и экспоненты.

http://en.wikipedia.org/wiki/IEEE_754-2008

Функция, имеющая float или double в числе аргументов, получает это значение через стек. Если функция возвращает float или double, она оставляет значение в регистре ST(0) - на вершине стека.

Простой пример

Рассмотрим простой пример:

double f (double a, double b)
{
	return a/3.14 + b*4.1;
};

Компилим с MSVC 2010:

CONST	SEGMENT
__real@4010666666666666 DQ 04010666666666666r	; 4.1
CONST	ENDS
CONST	SEGMENT
__real@40091eb851eb851f DQ 040091eb851eb851fr	; 3.14
CONST	ENDS
_TEXT	SEGMENT
_a$ = 8							; size = 8
_b$ = 16						; size = 8
_f	PROC
	push	ebp
	mov	ebp, esp
	fld	QWORD PTR _a$[ebp]
 
; состояние стека сейчас: ST(0) = _a
 
	fdiv	QWORD PTR __real@40091eb851eb851f
 
; состояние стека сейчас: ST(0) = результат деления _a на 3.13
 
	fld	QWORD PTR _b$[ebp]
 
; состояние стека сейчас: ST(0) = _b; ST(1) = результат деления _a на 3.13
 
	fmul	QWORD PTR __real@4010666666666666
 
; состояние стека сейчас: ST(0) = результат _b * 4.1; ST(1) = результат деления _a на 3.13
 
	faddp	ST(1), ST(0)
 
; состояние стека сейчас: ST(0) = результат сложения
 
	pop	ebp
	ret	0
_f	ENDP

FLD берет 8 байт из стека и загружает из в регистр ST(0), автоматически конвертируя во внутренний 80-битный формат.

FDIV делит содержимое ST(0) на число лежащее по адресу __real@40091eb851eb851f - там закодировано значение 3.14. Синтаксис ассемблера не поддерживает подобные числа, так что то что мы там видим, это шестандцатиричное представление числа 3.14.

После выполнения FDIV, в ST(0) остается результат деления.

Кстати, есть еще инструкция FDIVP, которая делит ST(1) на ST(0), выталкивает эти числа из стека и заталкивает результат. Если вы знаете язык Forth, то это как раз оно и есть - стековая машина.

Следующая FLD заталкивает в стек значение b.

После этого, в ST(1) перемещается результат деления, а в ST(0) теперь будет b.

Следующий FMUL умножает b из ST(0) на значение __real@4010666666666666 - там лежит число 4.1, и оставляет результат в ST(0).

Самая последняя инструкция FADDP складывает два значения из вершины стека, в ST(1) и затем выталкивает значение лежащее в ST(0), таким образом результат сложения остается на вершине стека в ST(0).

Функция должна вернуть результат в ST(0), так что больше ничего здесь не производится.

GCC 4.4.1 (с опцией -O3) генерит похожий код, хотя и с некоторой разницей:

                public f
f               proc near
 
arg_0           = qword ptr  8
arg_8           = qword ptr  10h
 
                push    ebp
                fld     ds:dbl_8048608 ; 3.14
 
; состояние стека сейчас: ST(0) = 3.13
 
                mov     ebp, esp
                fdivr   [ebp+arg_0]
 
; состояние стека сейчас: ST(0) = результат деления
 
                fld     ds:dbl_8048610 ; 4.1
 
; состояние стека сейчас: ST(0) = 4.1, ST(1) = результат деления
 
                fmul    [ebp+arg_8]
 
; состояние стека сейчас: ST(0) = результат умножения, ST(1) = результат деления
 
                pop     ebp
                faddp   st(1), st
 
; состояние стека сейчас: ST(0) = результат сложения
 
                retn
f               endp

Разница в том, что в стек сначала заталкивается 3.14 (в ST(0)), а затем значение из arg_0 делится на то что лежит в регистре ST(0).

FDIVR означает "Reverse Divide" - делить поменяв делитель и делимое местами. Точно такой же инструкции для умножения нет, потому что умножение - операция коммутативная, так что здесь остается FMUL.

FADDP не только складывает два значения, но также и выталкивает из стека одно значение. После этого, в ST(0) остается только результат сложения.

Этот кусок кода получен при помощи IDA, которая регистр ST(0) называет для краткости просто ST.

Передача чисел с плавающей запятой в аргументах

int main ()
{
	printf ("32.01 ^ 1.54 = %lf\n", pow (32.01,1.54));
 
	return 0;
}

pow() - это функция, которая возводит число в степень.

Посмотрим что у нас вышло (MSVC 2010):

CONST	SEGMENT
__real@40400147ae147ae1 DQ 040400147ae147ae1r	; 32.01
__real@3ff8a3d70a3d70a4 DQ 03ff8a3d70a3d70a4r	; 1.54
CONST	ENDS
 
_main	PROC
	push	ebp
	mov	ebp, esp
	sub	esp, 8  ; выделить место для первой переменной
	fld	QWORD PTR __real@3ff8a3d70a3d70a4
	fstp	QWORD PTR [esp]
	sub	esp, 8  ; выделить место для второй переменной
	fld	QWORD PTR __real@40400147ae147ae1
	fstp	QWORD PTR [esp]
	call	_pow
	add	esp, 8  ; "вернуть" место от одной переменной. 
 
; в локальном стеке сейчас зарезервировано 8 байт для нас. результат сейчас в ST(0)
 
	fstp	QWORD PTR [esp] ; перегрузить результат из ST(0) в локальный стек для printf()
	push	OFFSET $SG2651
	call	_printf
	add	esp, 12
	xor	eax, eax
	pop	ebp
	ret	0
_main	ENDP

FLD/FSTP перемешают переменные из/в сегмента данных в FPU-стек. pow() достает оба значения из FPU-стека и возвращает результат в ST(0). printf() берет 8 байт из стека и трактует как double.

Пример с сравнением

Попробуем теперь вот это:

double d_max (double a, double b)
{
	if (a>b)
		return a;
 
	return b;
};

Несмотря на кажущуюся простоту этой функции, понять как она работает будет чуть сложнее.

Вот что выдал MSVC 2010:

PUBLIC	_d_max
_TEXT	SEGMENT
_a$ = 8							; size = 8
_b$ = 16						; size = 8
_d_max	PROC
	push	ebp
	mov	ebp, esp
	fld	QWORD PTR _b$[ebp]
 
; состояние стека сейчас: ST(0) = _b
 
	fcomp	QWORD PTR _a$[ebp] ; compare _b (ST(0)) and _a and pop register
 
; стек теперь пустой
 
	fnstsw	ax
	test	ah, 5
	jp	SHORT $LN1@d_max
 
; we are here if a>b
 
	fld	QWORD PTR _a$[ebp]
	jmp	SHORT $LN2@d_max
$LN1@d_max:
	fld	QWORD PTR _b$[ebp]
$LN2@d_max:
	pop	ebp
	ret	0
_d_max	ENDP

Итак, FLD загружает _b в регистр ST(0).

FCOMP сравнивает содержимое ST(0) с тем что лежит в _a и выставляет биты C3/C2/C0 в регистре статуса FPU. Это 16-битный регистр отражающий текущее состояние FPU.

Итак, биты C3/C2/C0 выставлены, но, к сожалению, у процессора до Intel P6 нет инструкций условного перехода проверяющих эти биты. Возможно это исторически так сложилось (вспомните о том что FPU когда-то был вообще отдельным чипом). А у Intel P6 появились инструкции FCOMI/FCOMIP/FUCOMI/FUCOMIP - делающие тоже самое, только напрямую модифицирующие флаги ZF/PF/CF.

После этого, инструкция выдергивает одно значение из стека. Это отличает её от FCOM, которая просто сравнивает значения, оставляя стек в таком же состоянии.

FNSTSW копирует содержимое регистра статуса в AX. Биты C3/C2/C0 занимают позиции, соответственно, 14, 10, 8, в этих позициях они и остаются в регистре AX, и все они расположены в старшей части регистра - AH.

Если b>a в нашем случае, то биты C3/C2/C0 должны быть выставлены так: 0, 0, 0.

Если a>b, то биты будут выставлены: 0, 0, 1.

Если a=b, то биты будут выставлены так: 1, 0, 0.

После test ah, 5, бит C3 и C1 сбросится в ноль, на позициях 0 и 2 (внутри регистра AH) останутся соответственно C0 и C2.

Теперь немного о parity flag. Еще один замечательный рудимент эпохи.

http://en.wikipedia.org/wiki/Parity_flag

Там же в статье можно прочитать: One common reason to test the parity flag actually has nothing to do with parity. The FPU has four condition flags (C0 to C3), but they can not be tested directly, and must instead be first copied to the flags register. When this happens, C0 is placed in the carry flag, C2 in the parity flag and C3 in the zero flag. The C2 flag is set when e.g. incomparable floating point values (NaN or unsupported format) are compared with the FUCOM instructions.

Этот флаг выставляется в 1 если количество едениц в последнем результате - четно. И в ноль если - нечетно.

Таким образом, что мы имеем, флаг PF будет выставлен в еденицу, если C0 и C2 оба нули или оба еденицы. И тогда сработает последующий JP (jump if PF==1). Если мы вернемся чуть назад и посмотрим значения C3/C2/C0 для разных вариантов, то увидим, что условный переход JP сработает в двух случаях: если b>a или если a==b (ведь бит C3 уже "вылетел" после исполнения test ah, 5).

Дальше все просто. Если условный переход сработал, то FLD загрузит в ST(0) значение _b, а если не сработал, то загрузится _a и произойдет выход из функции.

Но это еще не все!

А теперь скомпилим все это в MSVC 2010 с опцией /Ox

_a$ = 8							; size = 8
_b$ = 16						; size = 8
_d_max	PROC
	fld	QWORD PTR _b$[esp-4]
	fld	QWORD PTR _a$[esp-4]
 
; состояне стека сейчас: ST(0) = _a, ST(1) = _b
 
	fcom	ST(1) ; compare _a and ST(1) = (_b)
	fnstsw	ax
	test	ah, 65					; 00000041H
	jne	SHORT $LN5@d_max
	fstp	ST(1)  ; copy ST(0) to ST(1) and pop register, leave (_a) on top
 
; состояние стека сейчас: ST(0) = _a
 
	ret	0
$LN5@d_max:
	fstp	ST(0)  ; copy ST(0) to ST(0) and pop register, leave (_b) on top
 
; состояние стека сейчас: ST(0) = _b
 
	ret	0
_d_max	ENDP

FCOM отличается от FCOMP тем что просто сравнивает значения и оставляет стек в том же состоянии. В отличие от предыдущего примера, операнды здесь в другом порядке. Поэтому и результат в C3/C2/C0 будет другим чем раньше:

Если a>b в нашем случае, то биты C3/C2/C0 должны быть выставлены так: 0, 0, 0.

Если b>a, то биты будут выставлены: 0, 0, 1.

Если a=b, то биты будут выставлены так: 1, 0, 0.

test ah, 65 как бы оставляет только два бита - C3 и C0. Они оба будут нулями, если a>b: в таком случае переход JNE не сработает. Далее имеется инструкция FSTP ST(1) - эта инструкция копирует значение ST(0) в указанный операнд и выдергивает одно значение из стека. В данном случае, она копирует ST(0) (где сейчас лежит _a) в ST(1). После этого на вершине стека два раза лежат _a. Затем одно значение выдергивается. После этого в ST(0) остается _a и функция завершается.

Условный переход JNE сработает в двух других случаях: если b>a или a==b. Скопируется ST(0) в ST(0), что как бы холостая операция, затем одно значение из стека вылетит и на вершине стека останется то что до этого лежало в ST(1) (то есть, _b). И функция завершится. Эта инструкция используется здесь видимо потому что в FPU нет инструкции которая просто выдергивает значение из стека и больше ничего.

Но и это еще не все.

GCC 4.4.1

d_max           proc near
 
b               = qword ptr -10h
a               = qword ptr -8
a_first_half    = dword ptr  8
a_second_half   = dword ptr  0Ch
b_first_half    = dword ptr  10h
b_second_half   = dword ptr  14h
 
                push    ebp
                mov     ebp, esp
                sub     esp, 10h
 
; переложим a и b в локальный стек:
 
                mov     eax, [ebp+a_first_half]
                mov     dword ptr [ebp+a], eax
                mov     eax, [ebp+a_second_half]
                mov     dword ptr [ebp+a+4], eax
                mov     eax, [ebp+b_first_half]
                mov     dword ptr [ebp+b], eax
                mov     eax, [ebp+b_second_half]
                mov     dword ptr [ebp+b+4], eax
 
; загружаем a и b в стек FPU
 
                fld     [ebp+a]
                fld     [ebp+b]
 
; текущее состояние стека: ST(0) - b; ST(1) - a
 
                fxch    st(1) ; эта инструкция меняет ST(1) и ST(0) местами
 
; текущее состояние стека: ST(0) - a; ST(1) - b
 
                fucompp            ; сравнить a и b и выдернуть из стека два значения, т.е., a и b
                fnstsw  ax         ; записать статус FPU в AX
                sahf               ; загрузить состояние флагов SF, ZF, AF, PF, и CF из AH
                setnbe  al                ; записать еденицу в AL если CF=0 и ZF=0
                test    al, al            ; AL==0 ?
                jz      short loc_8048453 ; да
                fld     [ebp+a]
                jmp     short locret_8048456
 
loc_8048453:
                fld     [ebp+b]
 
locret_8048456:
                leave
                retn
d_max           endp

FUCOMPP - это почти то же что и FCOM, только выкидывает из стека оба значения после сравнения, а также несколько иначе реагирует на "не числа".

Отступление:

http://en.wikipedia.org/wiki/NaN

http://ru.wikipedia.org/wiki/NaN

FPU умеет работать со специальными переменными, которые числами не являются. Это бесконечность, результат деления на ноль, и так далее. Нечисла бывают "тихие" и "сигнализирующие". С первыми можно продолжать работать и далее, а вот если вы попытаетесь совершить какую-то операцию с сигнализирующим нечислом, то сработает исключение. FCOM вызовет исключение если любой из операндов - нечисло. FUCOM же вызовет исключение только если один из операндов - сигнализирующее нечисло.

Далее мы видим SAHF - это довольно редкая инструкция в коде не использущим FPU. 8 бит из AH перекладываются в младшие 8 бит регистра статуса процессора в таком порядке: SF:ZF:?:AF:?:PF:?:CF <- AH.

Вспомним, что FNSTSW перегружает интересующие нас биты C3/C2/C0 в AH, и соответственно они будут в позициях 6, 2, 0 в регистре AH.

Иными словами, пара инструкций fnstsw ax / sahf перекладывает биты C3/C2/C0 в флаги ZF, PF, CF.

Теперь снова вспомним, какие значения бит C3/C2/C0 будут при каких результатах сравнения:

Если a больше b в нашем случае, то биты C3/C2/C0 должны быть выставлены так: 0, 0, 0.

Если a меньше b, то биты будут выставлены: 0, 0, 1.

Если a=b, то биты будут выставлены так: 1, 0, 0.

Иными словами, после инструкций FUCOMPP/FNSTSW/SAHF, мы получим такое состояние флагов:

Если a>b в нашем случае, то флаги будут выставлены так: ZF=0, PF=0, CF=0.

Если a меньше b, то флаги будут выставлены: ZF=0, PF=0, CF=1.

Если a равно b, то флаги будут выставлены так: ZF=1, PF=0, CF=0.

Инструкция SETNBE выставит в AL еденицу или ноль, в зависимости от флагов и условий. Это почти аналог JNBE, за тем лишь исключением, что SETcc (cc это condition code) выставляет 1 или 0 в AL, а Jcc делает переход или нет. SETNBE запишет еденицу если только CF=0 и ZF=0. Если это не так, то запишет 0.

CF будет 0 и ZF будет 0 одновременно только в одном случае, если a>b.

Тогда в AL будет записана еденица, послудющий условный переход JZ взят не будет, и функция вернет _a. В остальных случаях, функция вернет _b.

Но и это еще не конец.

GCC 4.4.1 с оптимизацией -O3

                public d_max
d_max           proc near
 
arg_0           = qword ptr  8
arg_8           = qword ptr  10h
 
                push    ebp
                mov     ebp, esp
                fld     [ebp+arg_0] ; _a
                fld     [ebp+arg_8] ; _b
 
; состояние стека сейчас: ST(0) = _b, ST(1) = _a
 
                fxch    st(1)
 
; состояние стека сейчас: ST(0) = _a, ST(1) = _b
 
                fucom   st(1) ; сравнить _a и _b
                fnstsw  ax
                sahf
                ja      short loc_8048448
                fstp    st ; записать ST(0) в ST(0) (холостая операция), выкинуть значение лежащее на вершине стека, оставить _b
                jmp     short loc_804844A
 
loc_8048448:
                fstp    st(1) ; записать _a в ST(0), выкинуть значение лежащее на вершине стека, оставить _a на вершине стека
 
loc_804844A:
                pop     ebp
                retn
d_max           endp

Почти все что здесь есть уже описано мною, кроме одного: использование JA после SAHF. Действительно, инструкции условных переходов "больше", "меньше", "равно" для сравнения беззнаковых чисел (JA, JAE, JB, JBE, JE/JZ, JNA, JNAE, JNB, JNBE, JNE/JNZ) проверяют только флаги CF и ZF. И биты C3/C2/C0 после сравнения перекладываются в эти флаги аккурат так, чтобы перечисленные инструкции переходов могли работать. JA сработает если CF и ZF обнулены.

Таким образом, перечисленные инструкции условного перехода можно использовать после инструкций FNSTSW/SAHF.

Вполне возможно что биты статуса FPU C3/C2/C0 преднамерено были размещены таким образом, чтобы переноситься на флаги процессора без перестановок.