Cтандарт RISC-V RV32I и математика с плавающей точкой
Никакое серьезное вычисление невозможно без чисел с плавающей запятой, но, как известно, стандарт RISC-V RV32I не содержит команд умножения и деления. Достаточно ли будет в софт-ядре реализовать стандарт RV32I, чтобы можно было вычислять что-то серьезное? В этой статье на примере RISC-V процессора YRV, описанного в книге Inside an Open-Source Processor, мы рассмотрим, как, используя компилятор GCC, рассчитать такие тригонометрические функции, как синус, косинус, тангенс, выведем на экран результат и даже нарисуем олдскульную синусоиду в VGA-режиме.
Софтовая математика
Когда смотришь на множество стандартов RISC-V, создается впечатление, что если ты не реализуешь в софт-ядре бо́льшую часть стандартов, компилятор ничего полезного для такого ядра собрать и не сможет. Поэтому попробуем скомпилировать следующий код для стандарта RV32I при помощи GCC и CLANG и посмотрим, что получится:
int main()
{
for (float i = 0; i < 10; i++)
{
for (float j = 0; j < 10; j++)
{
float k = i / j;
}
}
}В итоге команда
clang -target riscv32 -march=rv32i -S main.cвыдаст нам следующий результат (привожу его часть):
.LBB0_3: # Parent Loop BB0_1 Depth=1
# => This Inner Loop Header: Depth=2
lw a0, -20(s0)
lui a1, 266752
call __ltsf2
bgez a0, .LBB0_6
j .LBB0_4
.LBB0_4: # in Loop: Header=BB0_3 Depth=2
lw a0, -16(s0)
lw a1, -20(s0)
call __divsf3
sw a0, -24(s0)
j .LBB0_5
.LBB0_5: # in Loop: Header=BB0_3 Depth=2
lw a0, -20(s0)
lui a1, 260096
call __addsf3
sw a0, -20(s0)
j .LBB0_3
В случае сборки
gcc -c -march=rv32izbs -mabi=ilp32 -O0 -S main.s main.cрезультат будет следующим (также привожу часть):
.L4:
lw a1,-24(s0)
lw a0,-28(s0)
call __divsf3@plt
mv a5,a0
sw a5,-20(s0)
lla a5,.LC0
lw a1,0(a5)
lw a0,-24(s0)
call __addsf3@plt
mv a5,a0
sw a5,-24(s0)
.L3:
lla a5,.LC1
lw a1,0(a5)
lw a0,-24(s0)
call __ltsf2@plt
mv a5,a0
blt a5,zero,.L4
lla a5,.LC0
lw a1,0(a5)
lw a0,-28(s0)
call __addsf3@plt
mv a5,a0
sw a5,-28(s0)
И тот же код, но собранный для стандарта RV32G:
.L4:
flw fa4,-28(s0)
flw fa5,-24(s0)
fdiv.s fa5,fa4,fa5
fsw fa5,-20(s0)
flw fa4,-24(s0)
lla a5,.LC0
flw fa5,0(a5)
fadd.s fa5,fa4,fa5
fsw fa5,-24(s0)
.L3:
flw fa4,-24(s0)
lla a5,.LC1
flw fa5,0(a5)
flt.s a5,fa4,fa5
bne a5,zero,.L4
flw fa4,-28(s0)
lla a5,.LC0
flw fa5,0(a5)
fadd.s fa5,fa4,fa5
fsw fa5,-28(s0)
Нетрудно увидеть, что компилятор вместо нереализованной команды процессора вызывает соответствующие функции: для деления вместо исполнения команды fdiv.s вызывается функция divsf3, для сложения fadd.s — addsf3, сравнения flt.s — __ltsf2. В общем-то, логично. Можно заметить, что и регистры другие, но float входит в 32-битный РОН RV32I, поэтому в нашем случае это не принципиально. Для чистоты эксперимента повторим то же самое для типа int:
#RV32I
.L4:
lw a1,-24(s0)
lw a0,-28(s0)
call __divsi3@plt
#RV32G
.L4:
lw a4,-28(s0)
lw a5,-24(s0)
div a5,a4,a5
Жить можно: компилятор использует функции для команд, которых нет в стандарте, осталось найти их реализацию. Но есть один неприятный момент: на самом деле никакого межкомпиляторного стандарта соответствия ассемблерных инструкций и софтверных методов нет, если не считать стандартом сам GCC. В GCC эти методы определяются в библиотеке libgcc, а в Clang — в библиотеке builtins проекта compiler-rt, хотя имена методов совпадают.
Плавающая точка
И GCC и Clang поддерживают «софтовую» математику плавающей точки для RISC-V, да и, собственно, сам float прекрасно влезает в РОН RV32I. На мой взгляд, builtins Clang более удобен для изучения, так как большинство методов реализовано на C и лицензия более комфортна, а libgcc с первого взгляда непонятный, да и распространяется под лицензией GPLv3. Но мы вообще не будем пользоваться реализациями компиляторов, потому как в baremetal все должно быть под контролем.
В качестве библиотеки, реализующей математику с плавающей точкой, мы воспользуемся RVfplib — она проще для изучения и реализована на ассемблере, поэтому должна быть очень шустрой. Но есть одна неприятность: библиотека предназначена для архитектуры RV32IM, а у нас умножения нет.
Несложными экспериментами мы можем выяснить, что аппаратное умножение нужно нам лишь для одного метода — умножения. Заменим этот метод на любой из имеющихся у нас библиотек. Так как мы решили не пользоваться компиляторными библиотеками, возьмем библиотеку fp-soft, которая содержит C-реализации __mulsf3 (по-хорошему надо бы и лицензии на совместимость проверить). Небольшим недостатком RVfplib является то, что она не полностью поддерживает long double. Хоть на Марс мы не летим, это создает некоторые неудобства при использовании других библиотек.
fabs.o:fabs.c
$(GCC) $(CFLAGS) -o $@ $^
fp-precision.o:fp-precision.c
$(GCC) $(CFLAGS) -o $@ $^
fp-mul32.o:fp-mul32.c
$(GCC) $(CFLAGS) -o $@ $^
floatsisf.o:floatsisf.S
$(GCC) $(CFLAGS) -o $@ $^
addsf3.o:addsf3.S
$(GCC) $(CFLAGS) -o $@ $^
muldi3.o:muldi3.S
$(GCC) $(CCFLAGS) -o $@ $^
extendsfdf2.o:extendsfdf2.S
$(GCC) $(CFLAGS) -o $@ $^
fixsfsi.o:fixsfsi.S
$(GCC) $(CFLAGS) -o $@ $^
ltsf2.o:ltsf2.S
$(GCC) $(CFLAGS) -o $@ $^
ltdf2.o:ltdf2.S
$(GCC) $(CFLAGS) -o $@ $^
adddf3.o:adddf3.S
$(GCC) $(CFLAGS) -o $@ $^
divdf3.o:divdf3.S
$(GCC) $(CFLAGS) -o $@ $^
truncdfsf2.o:truncdfsf2.S
$(GCC) $(CFLAGS) -o $@ $^
gesf2.o: gesf2.S
$(GCC) $(CFLAGS) -o $@ $^
eqsf2.o:eqsf2.S
$(GCC) $(CFLAGS) -o $@ $^
eqdf2.o:eqdf2.S
$(GCC) $(CFLAGS) -o $@ $^
divsf3.o:divsf3.S
$(GCC) $(CFLAGS) -o $@ $^
floatsidf.o:floatsidf.S
$(GCC) $(CFLAGS) -o $@ $^
libfloat.a:muldi3.o floatsisf.o addsf3.o extendsfdf2.o fixsfsi.o ltsf2.o adddf3.o fp-mul32.o fp-precision.o divdf3.o truncdfsf2.o fabs.o gesf2.o eqsf2.o divsf3.o eqdf2.o ltdf2.o floatsidf.o
$(AR) $(ARFLAGS) -o $@ $^
В результате множества попыток сборки всего проекта получилась библиотечка, содержащая необходимые методы для работы с плавающей точкой. Все нужные ссылки на GitHub — в конце статьи.
Тригонометрия
Итак, если плавающая точка плывет, значит, она создает волну. Волны описываются функциями синуса или косинуса. А синус и косинус берется из таблиц Брадиса, но для baremetal это нереально. Кстати, в DOOM тригонометрия хоть и табличная, но для первичного расчета используется библиотечный синус, и в Wolfenstein 3D все аналогично.
for (i=0 ; i<5*FINEANGLES/4 ; i++)
{
// OPTIMIZE: mirror...
a = (i+0.5)*PI*2/FINEANGLES;
t = FRACUNIT*sin (a);
finesine[i] = t;
}
Возьмем CORDIC как алгоритм расчета синуса. Нам нужно проверить, как работает float, так что возьмем чужой алгоритм (Google говорит, что все используют этот Cordic). Сам алгоритм на целых числах, а значит, шустрый.
void cordic(int theta, int *s, int *c, int n)
{
int k, d, tx, ty, tz;
int x = cordic_1K, y = 0, z = theta;
n = (n > CORDIC_NTAB) ? CORDIC_NTAB : n;
for (k = 0; k < n; ++k)
{
d = z >> 31;
// get sign. for other architectures, you might want to use the more portable version
// d = z>=0 ? 0 : -1;
tx = x - (((y >> k) ^ d) - d);
ty = y + (((x >> k) ^ d) - d);
tz = z - ((cordic_ctab[k] ^ d) - d);
x = tx;
y = ty;
z = tz;
}
*c = x;
*s = y;
}
Поэтому подкинем граблей и добавим метод преобразования радиан в градусы. Все это будет люто неэффективно, но наша задача — проверить, что все будет работать. Да и какая плавающая точка без числа π.
void cordic_grad(int grad, int *sin, int *cos)
{
float p = grad * 3.14159265f / 180;
cordic((p * MUL), sin, cos, 32);
}
Вообще, я думал, что компилятор выбросит отсюда плавающую точку, так как MUL = 1073741824.000000f меньше MAX_INT, но, как мы видим ниже, компилятор использует «плавающее» деление (__divsf3) и умножение (__mulsf3).
cordic_grad:
.LFB1:
.cfi_startproc
addi sp,sp,-16
.cfi_def_cfa_offset 16
sw ra,12(sp)
sw s0,8(sp)
sw s1,4(sp)
.cfi_offset 1, -4
.cfi_offset 8, -8
.cfi_offset 9, -12
mv s0,a1
mv s1,a2
call __floatsisf@plt
lw a1,.LC0
call __mulsf3@plt
lw a1,.LC1
call __divsf3@plt
lw a1,.LC2
call __mulsf3@plt
call __fixsfsi@plt
li t1,652034048
lla t4,.LANCHOR0
Вывод на печать
Чтобы проверить правильность вычислений, хорошо бы увидеть результат, а для этого нам нужна какая-нибудь реализация printf. Для целых чисел и строк когда-то я использовал реализацию ee_printf из теста Coremark. Но та часть, которая предназначена для float, сложна в сборке (думаю, что ее никто для baremetal и не собирал). Здесь мне стало лень, и я дополнительно взял проект Nanoprintf, а точнее, метод npf_snprintf (…). Для вывода использовал уже имеющийся ee_printf (…), который уже отлажен на тесте Coremark.
char buf[32] = {0};
npf_snprintf(buf, sizeof(buf), " %d %0.6f", t, sin(t));
ee_printf("%s\n", buf);
Как я уже сказал, RVfplib не поддерживает long double, а Nanoprintf поддерживает, и это мешает линковке, поэтому эту поддержку я просто убрал.
#if NANOPRINTF_USE_FLOAT_FORMAT_SPECIFIERS == 1
// case 'L': out_spec->length_modifier = NPF_FMT_SPEC_LEN_MOD_LONG_DOUBLE; break;
#endif
Рисуем синус
Результирующий код выводит немного значений синусов, а также рисует два периода синусоиды. В реальном коде добавлены задержки для анимации.
void main()
{
clean(); //Чистим HEX
cls(); // Чистим экран
// Выводим значения синусов
for (int t = 0; t <= 360; t = t + 20)
{
char buf[32] = {0};
npf_snprintf(buf, sizeof(buf), " %d %0.6f", t, sin(t));
ee_printf("%s\n", buf);
}
// Рисуем два периода
for (int x = 0; x < 320; x++)
{
int grad = 720 * x / 320;
int y = (60.0f * sin(grad));
plot(x, 120 - y, 200);
}
wfi();
}
Железо
В качестве аппаратной платформы используем компьютер на базе RISC-V процессора YRV, описанный в статье «Ретро-компьютер уровня «Радио-86РК» с RISC-V процессором на плате OMDAZZ». Вся математика реализована программно, поэтому отладочной платы OMDAZZ недостаточно, на ней слишком мало памяти. Будем использовать плату DE0-CV, любезно предоставленную FPGA-Systems. Так как плата позволяет синтезировать больше 300 кБ оперативной памяти, сделаем графический режим, который реализован в виде VGA-like адаптера с разрешением 320×240, с доступом по адресу 0xA0000000. Все остальное аналогично варианту с OMDAZZ.
Селектор шины AHB-LITE для записи в видеопамять определяется константами VGA_BASE_х:
`define VGA_BASE_0 16'hA000
`define VGA_BASE_1 16'hA001 Пишем в видеопамять «если выполняется это условие»:
vga_wr_reg_0 <= mem_write && &mem_trans && (mem_addr[31:16] == `VGA_BASE_0);
vga_wr_reg_1 <= mem_write && &mem_trans && (mem_addr[31:16] == `VGA_BASE_1);
Основная проблема — это «узкая» память микросхемы Cyclone V, но зато ее много и она двухпортовая. Адрес пикселя вычисляется даже проще, чем в предыдущем проекте.
always @ (posedge clk)
begin
if(x== 799 && y==524)
pixel_addr<=0;
else
pixel_addr <= (((y>>1)<<8) + ((y>>1)<<6) + ((x+1)>>1));
end
X и y — координаты, получаемые от модуля vga.sv, стандартного для лабораторных работ Школы цифрового синтеза. Память двухпортовая, поэтому никаких проблем с чтением и записью нет.
Точка ставится следующим образом.
static char *VGA = (char *) 0xA0000000L;
void plot(int x, int y, char color) {
VGA[(y << 8) + (y << 6) + x] = color;
}
Работа с видеоадаптером подробно описана в книге Андре Ла Мот «Программирование игр для Microsoft Windows. Советы профессионала». Использование стандартного адреса режима 13h позволяет портировать все примеры из этой книги без особых проблем. Программа 8×8 Pixel ROM Font Editor позволяет сохранять шрифты в виде заголовочного файла font.h, который я разместил в библиотеке libyrv.a. Шрифт используется в методе ee_printf (…) от Coremark.
Результаты
Получилась у нас олдскульная ретросинусоида. Еще бы к ней трубку ЭЛТ…
Можно поиграть с точностью функции codric и посмотреть, как она влияет на вычисление. Для построения синусоиды достаточно n, равного 8 шагов, а 32 шага дают прекрасную точность.
А как же Soft FP от GCC?
Выполним следующий код (они никакого смысла не несет):
for (int t = 0; t <= 360; t = t + 20)
{
char buf[32] = {0};
float r = t/3;
float k = sin(t)*10;
npf_snprintf(buf, sizeof(buf), " %d %0.6f %0.6f ", t, sin(t),k+r);
ee_printf("%s\n", buf);
very_long_sleep();
}
И мы увидим, что дробная часть пропала:
Заменим функцию сложения из библиотеки libgcc. На самом деле, она просто собирается, нужно только добавить немного заголовочных файлов (double.h longlong.h op-1.h op-2.h op-4.h op-8.h op-common.h riscv-asm.h sfp-machine.h single.h soft-fp.h) и тогда наш Makefile прекрасно соберет эту функцию
addsf3.o:addsf3.c
$(GCC) $(CFLAGS) -o $@ $^Получим новую библиотеку, назовем ее gfloat (директория там же в репозитории). Сборка с ней выдаст необходимый нам результат:
Теперь все в порядке, но тогда почему бы не собрать библиотеку soft-fp из libgcc и не заморачиваться с какими-то внешними библиотеками? Пример библиотеки можно увидеть в директории gccfloat, но даже такой простой пример, как наш, собранный с этой библиотекой, занимает больше 32 кБ и не влезает в память нашего компьютера. Поэтому RVfplib важна даже не столько для скорости, сколько для размера. Сама библиотека soft-fp проверена (для чего размер памяти был увеличен до 128 кБ), но это уже совсем другая история.
Выводы
С помощью компилятора GCC стандарт RV32I позволяет реализовать вычисление математических функций с плавающей точкой одинарной точности. А алгоритмы реализации арифметических операций мозголомны и приятны глазу. На мой взгляд, на начальном этапе при изучении цифрового синтеза и стандарта RISC-V рациональнее сосредоточиться на базовом стандарте RV32I и больше уделять времени таким RISC-вещам, как конвейер, нежели делать однотактные процессоры, реализующие умножение и деление. Однако вопрос с лицензией на математические функции остается открытым. Как и обещал, прикладываю ссылку на репозиторий, с которым мы работали.
