2*2=4?

Меня давно мучал вопрос как объяснить gcc, что 2*2=4 - он весьма упёртый в этом деле, говорит мол 2*2=2, либо 4*4=4, а всё что кроме - мол не по стандарту. Где-то в инете видел упоминание о том, что научить его премудростям этой высшей математики таки можно, но как именно - тайна сакральная. По крайней мере мне лезть настолько глубоко в архитектуру этого монстра совершенно неохота. По-этому я пошёл другим путём - сделал ассемблерную вставку, реализующую то, что мне нужно. Сей код я тут оставлю - вдруг кому пригодится. Под катом реализация умножения int16_t*int16_t=int32_t, а заодно ещё квадратного корня uint16_t=sqrt(uint32_t) и корня суммы квадратов uint16_t=sqrt(int16_t^2 + int16_t^2).

Внимание! Я не гарантирую, что функции работают правильно на всём диапазоне. Если решите использовать эту фигню в своих проектах, обязательно проведите самостоятельное тестирование. Кроме того, нужно внимательно относиться к типу аргументов и результатов (знаковый/беззнаковый)
У умных дядек обычно принято в таких случаях указывать размер кода в словах и длительность выполнения в циклах, но я к ним не отношусь, по-этому могу только указать данные по исходным алгоритмам и время выполнения кода на осциллографу (atmega644PA на 1MHz).
Алгоритмы реализованы в виде макросов, если кому-то понадобятся функции - переделать вопрос пары минут.

1) Умножение int32_t a_result = int16_t a_x * int16_t a_y
Исходный алгоритм
Cycles : 17 + ret
Words : 13 + ret
Register usage: 11 registers
В моём варианте время выполнения 30мкс (встроенное умножение 32*32=32 - 66мкс).

#define MULS16X16_32(a_x, a_y, a_result)			\
{								\
	uint8_t tmp = 0;					\
	asm volatile (						\
	"clr	%[TMP]"			"\n\t" 			\
	"muls	%B[AX], %B[AY]"		"\n\t" 			\
	"movw	%C[RES], r0"		"\n\t" 			\
	"mul	%A[AX], %A[AY]"		"\n\t" 			\
	"movw	%A[RES], r0"		"\n\t" 			\
	"mulsu	%B[AX], %A[AY]"		"\n\t" 			\
	"sbc	%D[RES], %[TMP]"	"\n\t" 			\
	"add	%B[RES], r0"		"\n\t" 			\
	"adc	%C[RES], r1"		"\n\t" 			\
	"adc	%D[RES], %[TMP]"	"\n\t" 			\
	"mulsu	%B[AY], %A[AX]"		"\n\t" 			\
	"sbc	%D[RES], %[TMP]"	"\n\t" 			\
	"add	%B[RES], r0"		"\n\t" 			\
	"adc	%C[RES], r1"		"\n\t" 			\
	"adc	%D[RES], %[TMP]"	"\n\t" 			\
	"clr	r1"			"\n\t" 			\
		: [RES]"=&r" (a_result), [TMP]"=&r" (tmp)	\
		: [AX]"a" (a_x), [AY]"a" (a_y)			\
	);							\
}

2) Целочисленный квадратный корень uint16_t a_result = sqrt(uint32_t a_arg)
Исходный алгоритм
43 code words, 271 - 316 cycles
Register usage: 8 registers
В моём варианте время выполнения 280мкс (avr-libc-шный оптимизированный корень с плавающей точкой - 600-700мкс, из них примерно 500мкс, собственно, сам корень, а остальное - преобразование из целочисленных типов во float и обратно).
Внимание! Аргумент a_arg разрушается. Если он ещё нужен в дальнейшем - необходимо заводить временную переменную

#define SQRT32(a_arg, a_result)									\
{												\
	uint16_t tmp;										\
	asm volatile (										\
		"ldi   %B[TMP],0xc0"		"\n\t"						\
		"clr   %A[TMP]"			"; rotation mask in [TMP]\n\t"			\
		"ldi   %B[RES],0x40"		"\n\t"						\
		"sub   %A[RES],%A[RES]"		"; developing sqrt in [RES], C=0\n\t"		\
"_sq32_1_%=:"	"brcs  _sq32_2_%="		"; C --> Bit is always 1\n\t"			\
		"cp    %C[ARG],%A[RES]"		"\n\t"						\
		"cpc   %D[ARG],%B[RES]"		"; Does test value fit?\n\t"			\
		"brcs  _sq32_3_%="		"; C --> nope, bit is 0\n\t"			\
"_sq32_2_%=:"	"sub   %C[ARG],%A[RES]"		"\n\t"						\
		"sbc   %D[ARG],%B[RES]"		"; Adjust argument for next bit\n\t"		\
		"or    %A[RES],%A[TMP]"		"\n\t"						\
		"or    %B[RES],%B[TMP]"		"; Set bit to 1\n\t"				\
"_sq32_3_%=:"	"lsr   %B[TMP]"			"\n\t"						\
		"ror   %A[TMP]"			"; Shift right mask, C --> end loop\n\t"	\
		"eor   %B[RES],%B[TMP]"		"\n\t"						\
		"eor   %A[RES],%A[TMP]"		"; Shift right only test bit in result\n\t"	\
		"rol   %A[ARG]"			"; Bit 0 only set if end of loop\n\t"		\
		"rol   %B[ARG]"			"\n\t"						\
		"rol   %C[ARG]"			"\n\t"						\
		"rol   %D[ARG]"			"; Shift left remaining argument (C used at _sq32_1)\n\t" \
		"sbrs  %A[ARG],0"		"; Skip if 15 bits developed\n\t"		\
		"rjmp  _sq32_1_%="		"; Develop 15 bits of the sqrt\n\t"		\
		"brcs  _sq32_4_%="		"; C--> Last bits always 1\n\t"			\
		"cp    %A[RES],%C[ARG]"		"\n\t"						\
		"cpc   %B[RES],%D[ARG]"		"; Test for last bit 1\n\t"			\
		"brcc  _sq32_5_%="		"; NC --> bit is 0\n\t"				\
"_sq32_4_%=:"	"sbc   %B[ARG],%B[TMP]"		"; Subtract C (any value from 1 to 0x7f will do)\n\t" \
		"sbc   %C[ARG],%A[RES]"		"\n\t"						\
		"sbc   %D[ARG],%B[RES]"		"; Update argument for test\n\t"		\
		"inc   %A[RES]"			"; Last bit is 1\n\t"				\
"_sq32_5_%=:"	"lsl   %B[ARG]"			"; Only bit 7 matters\n\t"			\
		"rol   %C[ARG]"			"\n\t"						\
		"rol   %D[ARG]"			"; Remainder * 2 + C\n\t"			\
		"brcs  _sq32_6_%="		"; C --> Always round up\n\t"			\
		"cp    %A[RES],%C[ARG]"		"\n\t"						\
		"cpc   %B[RES],%D[ARG]"		"; C decides rounding\n\t"			\
"_sq32_6_%=:"	"adc   %A[RES],%B[TMP]"		"\n\t"						\
		"adc   %B[RES],%B[TMP]"		"; Round up if C (B[TMP]=0)\n\t"		\
			: [ARG]"=r" (a_arg), [RES]"=&d" (a_result), [TMP]"=&d" (tmp)		\
			: "[ARG]" (a_arg)							\
	);											\
}

3) Корень из суммы квадратов uint16_t a_result = sqrt((int16_t a_x)^2 + (int16_t a_x)^2)
Исходный алгоритм (примерно, без учёта подсчёта кормя)
; Size = 40 words
; Clock = 49 cycles
Register usage: 14 registers (в моём варианте, в исходном - хз)
В моём варианте время выполнения 320мкс.

#define SQRT_SQ_SUMM_16(a_x, a_y, a_result)							\
{												\
	uint16_t tmp16;										\
	uint32_t tmp32;										\
	asm volatile (										\
		"clr	%A[TMP16]"		"\n\t"						\
		"sbrs	%B[AX], 7"		"\n\t"						\
		"rjmp	_ss16_1_%="		"\n\t"						\
		"com	%A[AX]"			"\n\t"						\
		"com	%B[AX]"			"\n\t"						\
		"adc	%A[AX], %A[TMP16]"	"\n\t"						\
		"adc	%B[AX], %A[TMP16]"	"\n\t"						\
"_ss16_1_%=:"	"sbrs	%B[AY], 7"		"\n\t"						\
		"rjmp	_ss16_2_%="		"\n\t"						\
		"com	%A[AY]"			"\n\t"						\
		"com	%B[AY]"			"\n\t"						\
		"adc	%A[AY], %A[TMP16]"	"\n\t"						\
		"adc	%B[AY], %A[TMP16]"	"\n\t"						\
"_ss16_2_%=:"	"mul	%A[AY], %A[AY]"		"\n\t"						\
		"movw	%A[TMP32], r0"		"\n\t"						\
		"mul	%B[AY], %B[AY]"		"\n\t"						\
		"movw	%C[TMP32], r0"		"\n\t"						\
		"mul	%A[AY], %B[AY]"		"\n\t"						\
		"add	%B[TMP32], r0"		"\n\t"						\
		"adc	%C[TMP32], r1"		"\n\t"						\
		"adc	%D[TMP32], %A[TMP16]"	"\n\t"						\
		"add	%B[TMP32], r0"		"\n\t"						\
		"adc	%C[TMP32], r1"		"\n\t"						\
		"adc	%D[TMP32], %A[TMP16]"	"\n\t"						\
		"mul	%A[AX], %A[AX]"		"\n\t"						\
		"movw	%A[RES], r0"		"\n\t"						\
		"mul	%B[AX], %B[AX]"		"\n\t"						\
		"add	%A[TMP32], %A[RES]"	"\n\t"						\
		"adc	%B[TMP32], %B[RES]"	"\n\t"						\
		"adc	%C[TMP32], r0"		"\n\t"						\
		"adc	%D[TMP32], r1"		"\n\t"						\
		"mul	%A[AX], %B[AX]"		"\n\t"						\
		"add	%B[TMP32], r0"		"\n\t"						\
		"adc	%C[TMP32], r1"		"\n\t"						\
		"adc	%D[TMP32], %A[TMP16]"	"\n\t"						\
		"add	%B[TMP32], r0"		"\n\t"						\
		"adc	%C[TMP32], r1"		"\n\t"						\
		"adc	%D[TMP32], %A[TMP16]"	"\n\t"						\
		"clr	r1"			"\n\t"						\
		"ldi   %B[TMP16],0xc0"		"\n\t"						\
		"clr   %A[TMP16]"		"; rotation mask in [TMP16]\n\t"		\
		"ldi   %B[RES],0x40"		"\n\t"						\
		"sub   %A[RES],%A[RES]"		"; developing sqrt in [RES], C=0\n\t"		\
"_sq32_1_%=:"	"brcs  _sq32_2_%="		"; C --> Bit is always 1\n\t"			\
		"cp    %C[TMP32],%A[RES]"	"\n\t"						\
		"cpc   %D[TMP32],%B[RES]"	"; Does test value fit?\n\t"			\
		"brcs  _sq32_3_%="		"; C --> nope, bit is 0\n\t"			\
"_sq32_2_%=:"	"sub   %C[TMP32],%A[RES]"	"\n\t"						\
		"sbc   %D[TMP32],%B[RES]"	"; Adjust argument for next bit\n\t"		\
		"or    %A[RES],%A[TMP16]"	"\n\t"						\
		"or    %B[RES],%B[TMP16]"	"; Set bit to 1\n\t"				\
"_sq32_3_%=:"	"lsr   %B[TMP16]"		"\n\t"						\
		"ror   %A[TMP16]"		"; Shift right mask, C --> end loop\n\t"	\
		"eor   %B[RES],%B[TMP16]"	"\n\t"						\
		"eor   %A[RES],%A[TMP16]"	"; Shift right only test bit in result\n\t"	\
		"rol   %A[TMP32]"		"; Bit 0 only set if end of loop\n\t"		\
		"rol   %B[TMP32]"		"\n\t"						\
		"rol   %C[TMP32]"		"\n\t"						\
		"rol   %D[TMP32]"		"; Shift left remaining argument (C used at _sq32_1)\n\t" \
		"sbrs  %A[TMP32],0"		"; Skip if 15 bits developed\n\t"		\
		"rjmp  _sq32_1_%="		"; Develop 15 bits of the sqrt\n\t"		\
		"brcs  _sq32_4_%="		"; C--> Last bits always 1\n\t"			\
		"cp    %A[RES],%C[TMP32]"	"\n\t"						\
		"cpc   %B[RES],%D[TMP32]"	"; Test for last bit 1\n\t"			\
		"brcc  _sq32_5_%="		"; NC --> bit is 0\n\t"				\
"_sq32_4_%=:"	"sbc   %B[TMP32],%B[TMP16]"	"; Subtract C (any value from 1 to 0x7f will do)\n\t" \
		"sbc   %C[TMP32],%A[RES]"	"\n\t"						\
		"sbc   %D[TMP32],%B[RES]"	"; Update argument for test\n\t"		\
		"inc   %A[RES]"			"; Last bit is 1\n\t"				\
"_sq32_5_%=:"	"lsl   %B[TMP32]"		"; Only bit 7 matters\n\t"			\
		"rol   %C[TMP32]"		"\n\t"						\
		"rol   %D[TMP32]"		"; Remainder * 2 + C\n\t"			\
		"brcs  _sq32_6_%="		"; C --> Always round up\n\t"			\
		"cp    %A[RES],%C[TMP32]"	"\n\t"						\
		"cpc   %B[RES],%D[TMP32]"	"; C decides rounding\n\t"			\
"_sq32_6_%=:"	"adc   %A[RES],%B[TMP16]"	"\n\t"						\
		"adc   %B[RES],%B[TMP16]"	"; Round up if C (B[TMP16]=0)\n\t"		\
			: [RES]"=&d" (a_result), [TMP32]"=&r" (tmp32), [TMP16]"=&d" (tmp16)	\
			: [AX]"r" (a_x), [AY]"r" (a_y)						\
	);											\
}

Первоисточники алгоритмов:
1) AN201 (на сайте atmel`а его что-то нет, так что ссылка левая - http://faculty.capitol-college.edu/~andresho/tutor/Multimedia/AVR/HW_mult/avr201.htm)
2) http://members.chello.nl/j.beentjes3/Ruud/sqrt32avr.htm
3) http://elm-chan.org/docs/avrlib/sqrt32.S

UPD (02.02.2001): Пофиксил ошибки в типах аргументов, плюс улучшение читабельности и правка комментариев, где они есть.
























Смотрите также:

Вам это будет интересно!

  1. О работе с массивами