Мучаем 5576ХС4Т - часть 'h2E - делим отрезок на равные части
Часть 0 - покупаем, паяем, ставим драйвера и софт
Часть 1 - что это вообще за зверь?
Часть 2 - наша первая схема!
Часть 3 - кнопочки и лампочки
Часть 4 - делитель частоты
Часть 5 - подавление дребезга кнопки
Часть 6 - заканчиваем кнопочки и лампочки
Часть 7 - счетчики и жаба
Часть 8 - передатчик UART
Часть 9 - Hello, wolf!
Часть 'hA - приёмник UART
Часть 'hB - UART и жаба
Часть 'hC - полудуплексный UART.
Часть 'hD - МКО (МКИО, Mil-Std 1553) для бедных, введение.
Часть 'hE - приёмопередатчик МКО "из подручных материалов" (в процессе)
Часть 'hF - модуль передатчика МКО
Часть 'h10 - передатчик сообщений МКО
Часть 'h20 - работа с АЦП ADC124s051
Часть 'h21 - преобразование двоичного кода в двоично-десятичный (BCD)
Часть 'h22 - Bin2Bcd с последовательной выдачей данных
Часть 'h23 - перемножитель беззнаковых чисел с округлением
Часть 'h24 - перемножитель беззнаковых чисел, реализация
Часть 'h25 - передаём показания АЦП на компьютер
Часть 'h26 - работа над ошибками (быстрый UART)
Часть 'h27 - PNG и коды коррекции ошибок CRC32
Часть 'h28 - передатчик изображения PNG
Часть 'h29 - принимаем с ПЛИС изображение PNG
Часть 'h2A - ZLIB и коды коррекции ошибок Adler32
Часть 'h2B - ускоряем Adler32
Часть 'h2C - формирователь потока Zlib
Часть 'h2D - передаём сгенерированное PNG-изображение
Часть 'h2E - делим отрезок на равные части
Часть 'h2F - знаковые умножители, тысячи их!
Часть 'h30 - вычислитель множества Мандельброта
Часть 'h31 - ускоренные сумматоры
Часть 'h32 - ускоренные счётчики (делаем часы)
Часть 'h33 - ускоряем ВСЁ
Часть 'h34 - ускоренные перемножители
Часть 'h35 - умножители совсем просто
Часть 'h36 - уравновешенный четверичный умножитель
Мы хотим построить на ПЛИС множество Мандельброта - мне кажется, это хорошая тренировка как работы с большими (1920х1080) изображениями, так и выполнения арифметических операций на ПЛИС, включая работу с комплексными числами, хотя "в явном виде" вводить их здесь необязательно.
И действительную, и мнимую часть комплексного числа мы представим 16-битными знаковыми числами с фиксированной точкой с 2 битами перед запятой. Иногда такой формат называют Q2.14 (Q намекает, что мы работаем с рациональными числами, 2 бита перед запятой, 14 бит - после неё). В таком формате можно представить числа в диапазоне от -2 до 2-2-14. Для построения множества Мандельброта этого диапазона вполне достаточно - если число выходит за его границы - назад оно вернуться уже не сможет (последовательность получается расходящейся), и множеству Мандельброта данная точка не принадлежит.
Мы хотим, чтобы по оси Y помещались значения от -1 до +1, а по оси X - от -2 до ≈1,56 (для соблюдения пропорций).
И первой же задачей становится вычисление выражений
y0 = (j-540.0) / 540.0;
x0 = -2.0 + i / 540.0;
то есть, пересчитать координаты очередного пикселя (i, j) в комплексное число z0=x0 + iy0, которое этому пикселю соответствует.
Оказывается, что если мы идём от пикселя к пикселю слева направо, сверху вниз, то на каждом шаге нам достаточно двух сложений и проверки знака.
К сожалению, не знаю, есть ли специальные обозначения, чтобы не путаться, когда мы думаем о 16-битном числе как принимающем значения от -2 до 2-2-14, а когда - в качестве целого числа в диапазоне от -32768 до +32767. В формулах выше мы имели в виду первое, т.е j = 0 даст y0 = -1, а j = 1080 приведёт к y0 = +1.
Но если записать x0, y0 как целочисленные, то формулы надо домножить на 214 = 16384:
y0 = j * 16384 / 540 - 16384;
x0 = -32768 + i * 16384 / 540;
Немножко сократим дроби:
y0 = j * 4096 / 135 - 16384;
x0 = -32768 + i * 4096 / 135;
В большинстве языков программирования, целочисленное деление "округляет вниз", т.е 8 / 9 даст 0. Чтобы заставить его "округлять до ближайшего целого", надо к числителю прибавить половинку знаменателя:
y0 = (j * 4096 + 67) / 135 - 16384;
x0 = -32768 + (i * 4096 + 67) / 135;
(кстати, деление - очень мерзопакостная операция ещё и в том, что округляет скорее не "вниз", а "к нулю", т.е в большинстве реализаций -8 / 9 также даст 0, что нарушает "изотропность" операций - вокруг нуля имеем "зону нечувствительности" расширенную вдвое. Здесь мы по счастью делим только положительные числа - одной проблемой меньше...)
Для осуществления на 32-разрядном процессоре, у которого есть команда деления - это уже неплохие формулы. Они дадут правильный результат, в пределах тех 16 бит, в которых мы храним результат.
Но когда у нас ПЛИС, то даже такие формулы кажутся слишком сложными: нам нужно делать умножение с сохранением результата в 32-битном аккумуляторе, затем ещё и провести деление - мягко говоря не самая любимая операция для реализации "в железе" (причём как ни странно, целочисленное деление оказывается в чём-то даже сложнее деления с плавающей запятой, поскольку там оба операнда фактически принимают значения от 1 до 2, а тут - какие угодно!).
По счастью, наша задача чуточку проще: мы знаем, что первоначально каждая из переменных i,j равна нулю, а затем мы прибавляем единичку к i во время движения по строке, после чего снова зануляем её и прибавляем единичку к j при переходе на новую строку.
Так что всё, что нам нужно - это правильно обновлять значения x0, y0 на каждом шаге.
Немножко перепишем выражения:
y0 = 30 * j + (46 * j + 67) / 135 - 16384;
x0 = 30 * i + (46 * j + 67) / 135 - 32768;
Мы видим, что на каждом шаге к переменной прибавляется 30 (за счёт первого слагаемого), и ещё должно добавляться 46 / 135 ≈ 0,34, но поскольку мы считаем в целых числах, то эти 0,34 будут "пропадать", но достаточно накопившись, дадут прибавление ещё единицы.
Собственно, если честно протабулировать значения y0, мы получаем:
-16384
-16354 (+30)
-16323 (+31)
-16293 (+30)
-16263 (+30)
-16232 (+31)
Получается, что на каждом шаге нам нужно лишь "принять решение" - прибавлять ли единицу? Само прибавление осуществить легко - в сумматоре есть вход переноса с младших разрядов - а ведь именно это у нас и происходит по сути :)
Чтобы принять решение, мы по сути храним "остаток от деления". Во время инициализации у нас уже осталась дробь 67/135, которую мы округлили до нуля. Знаменатель всегда будет 135 - его мы хранить не должны, он "сам собой подразумевается". А вот значение 67 мы запомним.
Пошёл первый шаг, и теперь дробь стала (67 + 46) / 135 = 113 / 135. Числитель всё ещё меньше знаменателя, поэтому единицу не прибавляем. Запоминаем значение 113.
На следующем шаге мы снова прибавляем к числителю 46, получая дробь (113 + 46) / 135 = 159 / 135. Числитель стал больше знаменателя - значит, мы прибавляем единицу к значению y0, а из нашего "остатка" эту единицу надо убрать, для чего мы вычитаем из числителя знаменатель, т.е получаем (159 - 135) / 135 = 24 / 135.
Вот так и работаем - ничего кроме сложений, вычитаний (а это одно и то же по сути) и сравнений нам не нужно. Есть целое семейство методов, схожих с этим, с помощью которых рисуют прямые, окружности, кривые Безье и прочие "графические примитивы" на растровом изображении.
Осталось только максимально "приспособить" эту задачу под специфику ПЛИС. Сравнение двух чисел "в лоб" - операция довольно неприятная, лучше было бы, если достаточно было посмотреть на знак числа. Для этого нужно "сдвинуть" дробь, чтобы именно "пересечение нуля" означало необходимость добавить единицу.
То есть, работа должна начинаться не с 67 / 135, а с 67 / 135 - 1 = -68 / 135. Прибавляем 46 первый раз - получаем -22 / 135 (пока ещё минус), прибавляем 46 второй раз - получаем 24 / 135 - это положительное число, поэтому прибавляем единицу к y0, вычитаем из числителя 135 и приходим к -111 / 135.
Вероятно, лучше всего будет ещё и "инвертировать" нашу дробь, чтобы старший разряд нашего "остатка" (в дополнительном коде появление единицы в старшем разряде означает - число стало отрицательным) напрямую подключался к входу переносу сумматора, не требуя дополнительной логики. То есть, начинаем мы с 68 / 135, на каждом шагу вычитаем 46, а если результат стал отрицательным - прибавляем единицу к y0, а к числителю прибавляем 135.
Наконец, замечаем, что не нужно спешить с прибавлением 135. Пусть по окончании шага результат остаётся отрицательным - просто в следующий раз, если числитель отрицательный, то вместо вычитания 46 мы будем прибавлять к нему 89. Это, к тому же, уменьшит диапазон чисел, который нам надо представлять в этом "остатке". Теперь самым большим числом может оказаться только 88 (ведь мы прибавляем значения лишь к отрицательным числам), а самым маленьким - число (-46), поскольку в худшем случае мы однажды придём к нулю, и в следующий раз вычтем из него 46. В нашем примере хватит 8-битного числа со знаком.
Пора приступать к кодированию на Verilog!
Первым делом, поиграемся с функциями и выражениями, исполняющимися на этапе компиляции:
module EqualPartsDivider (input clk, input sclr, input ce, output [15:0] Q);
parameter InitValue = -16384; //-16384
parameter Numerator = 32768;
parameter Denominator = 1080;
localparam CommonDivider = GCD (Numerator, Denominator);
localparam SimpleNum = Numerator / CommonDivider;
localparam SimpleDenom = Denominator / CommonDivider;
localparam WholePart = SimpleNum / SimpleDenom;
localparam Residue = SimpleNum - WholePart * SimpleDenom;
localparam RInitValue = (SimpleDenom + 1) / 2 - Residue;
localparam DenomMinusRes = SimpleDenom - Residue;
function integer GCD (input integer A, input integer B);
integer r;
begin
while (B != 0) begin
r = A % B;
A = B;
B = r;
end
GCD = A;
end
endfunction
initial begin
$display ("CommonDivider=%d", CommonDivider);
$display ("SimpleNum=%d", SimpleNum);
$display ("SimpleDenom=%d", SimpleDenom);
$display ("WholePart=%d", WholePart);
$display ("Residue=%d", Residue);
$display ("RInitValue=%d", RInitValue);
$display ("DenomMinusRes=%d", DenomMinusRes);
end
endmodule
Мы решили "довериться" компилятору в плане нахождения всех значений. Мы задаём ему 3 параметра: начальное значение, числитель и знаменатель, т.е говорим, что хотим получать значения выражения
Q = InitValue + Round (Numenator * i / Denominator);
для последовательных значений i.
Дальше это дело компилятора сообразить, где что надо прибавлять или вычитать.
Мы в кои-то веки написали функцию нахождения общего делителя (алгоритм Эвклида), и применять её мы будем только на этапе компиляции.
(вообще, функции могут быть доступны и в рантайме, но с этим связано множество проблем - иногда код просто не будет синтезироваться, поскольку не знает, сколько итераций будет внутри, а иногда даже синтезируется, но окажется чрезвычайно громоздким и медленным)
Затем мы сократили числитель и знаменатель, обозвав эти значения SimpleNum и SimpleDenom.
После этого мы выделили целую часть WholePart, и остаток Residue, и нашли исходное значение для нашего вспомогательного регистра, RInitValue (об этом чуточку позже). И напоследок, находим, что мы должны прибавить к остатку, когда он с прошлого раза остался отрицательным, это число мы назвали DenomMinusRes (знаменатель минус остаток).
Чтобы проверить правильность данных операций, мы решили задействовать несинтезируемую конструкцию $display - Квартус её тупо игнорирует, но есть симуляторы verilog, которые воспримут её и "напишут в консоль" то, что мы попросим. В том числе, есть онлайн-ресурсы, например, https://www.tutorialspoint.com/compile_verilog_online.php
(многие ютуберы показывают примеры на EDAboard, но там требуется регистрация, причём они ненавидят бесплатные почтовые сервисы и почтовые ящики, лежащие на серверах "конкурентов". Ну их нафиг...)
Вот что мы получаем:
$iverilog -o main *.v
$vvp main
CommonDivider= 8
SimpleNum= 4096
SimpleDenom= 135
WholePart= 30
Residue= 46
RInitValue= 22
DenomMinusRes= 89
135, 30, 46, 89 - эти числа мы совсем недавно использовали :)
Разберём заголовок модуля:
module EqualPartsDivider (input clk, input sclr, input ce, output reg [15:0] Q = 1'b0);
Интерфейс довольно простой: у нас есть вход тактовой частоты clk, вход для синхронного сброса sclr и вход разрешения счёта ce. Результирующее число мы получаем на выходе Q, причём оно устанавливается по положительному фронту clk во время подачи сигнала разрешения счёта ce, и сохраняет верное значение до последующих действий (сброса или счёта). Таким образом, "отлавливать", когда значение на выходе верное, нам не нужно.
Выход мы сразу определили как регистр. Ещё один регистр нужен для хранения "остатка". Его ширина может варьироваться в зависимости от дроби, которую мы пытаемся посчитать. Максимальное отрицательное значение не может превысить числа (-Residue), максимальное положительное значение - числа DenomMinusRes. Мы можем найти требуемую ширину регистра следующим образом:
localparam ResidueWidth1 = `CLOG2((DenomMinusRes + 1));
localparam ResidueWidth2 = `CLOG2(Residue);
localparam ResidueWidth = 1 + ((ResidueWidth1 > ResidueWidth2)? ResidueWidth1 : ResidueWidth2);
Макрос `CLOG2 мы сделали в части 6 - чрезвычайно полезная штука! Но буквально сейчас на ней обжёгся - поначалу у меня были строки
localparam ResidueWidth1 = `CLOG2(DenomMinusRes + 1) + 1;
localparam ResidueWidth2 = `CLOG2(Residue) + 1;
localparam ResidueWidth = (ResidueWidth1 > ResidueWidth2)? ResidueWidth1 : ResidueWidth2;
и "+1" в конце каждой строки у меня "сожралось", потому что `CLOG2 - это прямая подстановка выражения
(x <= 2) ? 1 : \
(x <= 4) ? 2 : \
(x <= 8) ? 3 : \
(x <= 16) ? 4 : \
(x <= 32) ? 5 : \
(x <= 64) ? 6 : \
(x <= 128)? 7 : \
(x <= 256)? 8 : \
(x <= 512)? 9 : \
(x <= 1024)? 10 : \
(x <= 2048)? 11 : \
(x <= 4096)? 12 : \
(x <= 8192)? 13 : \
(x <= 16384)? 14 : \
(x <= 32768)? 15 : \
(x <= 65536)? 16 : \
(x <= 'h20000)? 17 : \
(x <= 'h40000)? 18 : \
(x <= 'h80000)? 19 : \
(x <= 'h100000)? 20 : \
(x <= 'h200000)? 21 : \
(x <= 'h400000)? 22 : \
(x <= 'h800000)? 23 : \
(x <= 'h1000000)? 24 : \
(x <= 'h2000000)? 25 : \
(x <= 'h4000000)? 26 : \
(x <= 'h8000000)? 27 : \
(x <= 'h10000000)? 28 : \
(x <= 'h20000000)? 29 : \
(x <= 'h40000000)? 30 : \
(x <= 'h80000000)? 31 : \
(x <= 'h0FFFFFFFF)? 32 : \
-1
в строку. Единица прибавилась к последнему варианту. Чтобы такого не происходило, надо скобки ставить, чем больше - тем лучше! Всё как в препроцессоре C.
В нашем примере, ширина регистра составит 8 бит. Мы будем рассматривать его как число в дополнительном коде, т.е значения от -128 до +127.
Напишем, наконец, логику нашего модуля! На "чистом верилоге" она очень проста:
reg [ResidueWidth - 1 : 0] D = 1'b0; //Decision
wire IsNeg = D[ResidueWidth - 1];
assign DD = D;
always @(posedge clk) if (ce | sclr) begin
D <= sclr? RInitValue : IsNeg? D + DenomMinusRes : D - Residue;
Q <= sclr? InitValue : Q + WholePart + IsNeg;
end
И действительно, этот код заработает так, как нам надо. Но результат синтеза меня не радует: для наших параметров получилось 78 ЛЭ. Точнее, само по себе число 78 - не такое уж страшное (если только это не заголовок Zlib), хуже другое: в логах мы видим, что в этих двух строках компилятор усмотрел аж 4 сумматора!
Как всегда, он не понял, что прибавление единички не требует дополнительного сумматора, ведь у нас есть вход переноса, ну ё-моё! И не понял, что можно применить один сумматор в первой строке, просто подставить ему разные значения.
Так что нам ничего не остаётся, как вставить сумматоры ручками:
reg [ResidueWidth - 1 : 0] D = 1'b0; //Decision
assign DD = D;
wire IsNeg = D[ResidueWidth - 1];
wire [ResidueWidth - 1 : 0] ValueToAdd = IsNeg? DenomMinusRes : (-Residue);
wire [ResidueWidth - 1: 0 ] DResult;
lpm_add_sub ResidueAdder ( .dataa (D),
.datab (ValueToAdd),
.result (DResult) );
defparam
ResidueAdder.lpm_direction = "ADD",
ResidueAdder.lpm_hint = "ONE_INPUT_IS_CONSTANT=NO,CIN_USED=NO",
ResidueAdder.lpm_representation = "SIGNED",
ResidueAdder.lpm_type = "LPM_ADD_SUB",
ResidueAdder.lpm_width = ResidueWidth;
wire [15:0] QResult;
lpm_add_sub MainAdder ( .dataa (Q),
.datab (WholePart),
.cin (IsNeg),
.result (QResult) );
defparam
MainAdder.lpm_direction = "ADD",
MainAdder.lpm_hint = "ONE_INPUT_IS_CONSTANT=YES,CIN_USED=YES",
MainAdder.lpm_representation = "SIGNED",
MainAdder.lpm_type = "LPM_ADD_SUB",
MainAdder.lpm_width = 16;
always @(posedge clk) if (ce | sclr) begin
D <= sclr? RInitValue : DResult;
Q <= sclr? InitValue : QResult;
end
Теперь мы получаем те же самые результаты, но используя лишь 50 ЛЭ, а что даже важнее - мы вдвое сократили комбинаторные пути, так что этот модуль спокойненько компилится и на медленном кристалле "-3", при заданной частоте в 80 МГц (даже остался запас по частоте аж в 10 МГц). Исходная схема уже ругалась...
Под спойлером - код модуля целиком (включая заголовок и толпу локальных параметров):
[Spoiler (click to open)]
`include "math.v"
module EqualPartsDivider (input clk, input sclr, input ce, output reg [15:0] Q = 1'b0, output [7:0] DD);
parameter InitValue = -16384; //-16384
parameter Numerator = 32768;
parameter Denominator = 1080;
localparam CommonDivider = GCD (Numerator, Denominator);
localparam SimpleNum = Numerator / CommonDivider;
localparam SimpleDenom = Denominator / CommonDivider;
localparam WholePart = SimpleNum / SimpleDenom;
localparam Residue = SimpleNum - WholePart * SimpleDenom;
localparam RInitValue = (SimpleDenom + 1) / 2 - Residue;
localparam DenomMinusRes = SimpleDenom - Residue;
localparam ResidueWidth1 = `CLOG2(DenomMinusRes + 1);
localparam ResidueWidth2 = `CLOG2(Residue);
localparam ResidueWidth = 1 + ((ResidueWidth1 > ResidueWidth2)? ResidueWidth1 : ResidueWidth2);
reg [ResidueWidth - 1 : 0] D = 1'b0; //Decision
assign DD = D;
wire IsNeg = D[ResidueWidth - 1];
wire [ResidueWidth - 1 : 0] ValueToAdd = IsNeg? DenomMinusRes : (-Residue);
wire [ResidueWidth - 1: 0 ] DResult;
lpm_add_sub ResidueAdder ( .dataa (D),
.datab (ValueToAdd),
.result (DResult) );
defparam
ResidueAdder.lpm_direction = "ADD",
ResidueAdder.lpm_hint = "ONE_INPUT_IS_CONSTANT=NO,CIN_USED=NO",
ResidueAdder.lpm_representation = "SIGNED",
ResidueAdder.lpm_type = "LPM_ADD_SUB",
ResidueAdder.lpm_width = ResidueWidth;
wire [15:0] QResult;
lpm_add_sub MainAdder ( .dataa (Q),
.datab (WholePart),
.cin (IsNeg),
.result (QResult) );
defparam
MainAdder.lpm_direction = "ADD",
MainAdder.lpm_hint = "ONE_INPUT_IS_CONSTANT=YES,CIN_USED=YES",
MainAdder.lpm_representation = "SIGNED",
MainAdder.lpm_type = "LPM_ADD_SUB",
MainAdder.lpm_width = 16;
always @(posedge clk) if (ce | sclr) begin
D <= sclr? RInitValue : DResult;
Q <= sclr? InitValue : QResult;
end
function integer GCD (input integer A, input integer B);
integer r;
begin
while (B != 0) begin
r = A % B;
A = B;
B = r;
end
GCD = A;
end
endfunction
initial begin
$display ("CommonDivider=%d", CommonDivider);
$display ("SimpleNum=%d", SimpleNum);
$display ("SimpleDenom=%d", SimpleDenom);
$display ("WholePart=%d", WholePart);
$display ("Residue=%d", Residue);
$display ("RInitValue=%d", RInitValue);
$display ("DenomMinusRes=%d", DenomMinusRes);
$display ("ResidueWidth=%d", ResidueWidth);
end
endmodule
Ну и объясним напоследок значение RInitValue - оно сюда внесено "с опережением на шаг"! Т.е в наших размышлениях мы начинали со значения 68 (половинка от знаменателя, округлённая "вверх"), а здесь сразу же вычли ещё 46, получив значение 22. Это позволяет сумматору для Q брать значение не с выхода сумматора для D (как мы делали в своих устных вычислениях - сначала вычитаем здесь, потом, если получили отрицательное значение - прибавим здесь единицу), а с выхода регистра. А для этого нужно, чтобы к началу такта в регистре уже хранилось правильное значение. По счастью, в этой задаче мы заранее знали, как его посчитать - вот и посчитали заблаговременно!
Для нашего генератора множества Мандельброта нужно два таких модуля, с разными параметрами. Заметим, что они не только вычисляют выражение, но и служат регистрами для хранения текущего Z0.
Смотрите в следующей части - "аппаратная" реализация выражения Z = Z2 + Z0, с проверкой на переполнение.
Poll Целочисленная арифметика
Часть 1 - что это вообще за зверь?
Часть 2 - наша первая схема!
Часть 3 - кнопочки и лампочки
Часть 4 - делитель частоты
Часть 5 - подавление дребезга кнопки
Часть 6 - заканчиваем кнопочки и лампочки
Часть 7 - счетчики и жаба
Часть 8 - передатчик UART
Часть 9 - Hello, wolf!
Часть 'hA - приёмник UART
Часть 'hB - UART и жаба
Часть 'hC - полудуплексный UART.
Часть 'hD - МКО (МКИО, Mil-Std 1553) для бедных, введение.
Часть 'hE - приёмопередатчик МКО "из подручных материалов" (в процессе)
Часть 'hF - модуль передатчика МКО
Часть 'h10 - передатчик сообщений МКО
Часть 'h20 - работа с АЦП ADC124s051
Часть 'h21 - преобразование двоичного кода в двоично-десятичный (BCD)
Часть 'h22 - Bin2Bcd с последовательной выдачей данных
Часть 'h23 - перемножитель беззнаковых чисел с округлением
Часть 'h24 - перемножитель беззнаковых чисел, реализация
Часть 'h25 - передаём показания АЦП на компьютер
Часть 'h26 - работа над ошибками (быстрый UART)
Часть 'h27 - PNG и коды коррекции ошибок CRC32
Часть 'h28 - передатчик изображения PNG
Часть 'h29 - принимаем с ПЛИС изображение PNG
Часть 'h2A - ZLIB и коды коррекции ошибок Adler32
Часть 'h2B - ускоряем Adler32
Часть 'h2C - формирователь потока Zlib
Часть 'h2D - передаём сгенерированное PNG-изображение
Часть 'h2E - делим отрезок на равные части
Часть 'h2F - знаковые умножители, тысячи их!
Часть 'h30 - вычислитель множества Мандельброта
Часть 'h31 - ускоренные сумматоры
Часть 'h32 - ускоренные счётчики (делаем часы)
Часть 'h33 - ускоряем ВСЁ
Часть 'h34 - ускоренные перемножители
Часть 'h35 - умножители совсем просто
Часть 'h36 - уравновешенный четверичный умножитель
Мы хотим построить на ПЛИС множество Мандельброта - мне кажется, это хорошая тренировка как работы с большими (1920х1080) изображениями, так и выполнения арифметических операций на ПЛИС, включая работу с комплексными числами, хотя "в явном виде" вводить их здесь необязательно.
И действительную, и мнимую часть комплексного числа мы представим 16-битными знаковыми числами с фиксированной точкой с 2 битами перед запятой. Иногда такой формат называют Q2.14 (Q намекает, что мы работаем с рациональными числами, 2 бита перед запятой, 14 бит - после неё). В таком формате можно представить числа в диапазоне от -2 до 2-2-14. Для построения множества Мандельброта этого диапазона вполне достаточно - если число выходит за его границы - назад оно вернуться уже не сможет (последовательность получается расходящейся), и множеству Мандельброта данная точка не принадлежит.
Мы хотим, чтобы по оси Y помещались значения от -1 до +1, а по оси X - от -2 до ≈1,56 (для соблюдения пропорций).
И первой же задачей становится вычисление выражений
y0 = (j-540.0) / 540.0;
x0 = -2.0 + i / 540.0;
то есть, пересчитать координаты очередного пикселя (i, j) в комплексное число z0=x0 + iy0, которое этому пикселю соответствует.
Оказывается, что если мы идём от пикселя к пикселю слева направо, сверху вниз, то на каждом шаге нам достаточно двух сложений и проверки знака.
К сожалению, не знаю, есть ли специальные обозначения, чтобы не путаться, когда мы думаем о 16-битном числе как принимающем значения от -2 до 2-2-14, а когда - в качестве целого числа в диапазоне от -32768 до +32767. В формулах выше мы имели в виду первое, т.е j = 0 даст y0 = -1, а j = 1080 приведёт к y0 = +1.
Но если записать x0, y0 как целочисленные, то формулы надо домножить на 214 = 16384:
y0 = j * 16384 / 540 - 16384;
x0 = -32768 + i * 16384 / 540;
Немножко сократим дроби:
y0 = j * 4096 / 135 - 16384;
x0 = -32768 + i * 4096 / 135;
В большинстве языков программирования, целочисленное деление "округляет вниз", т.е 8 / 9 даст 0. Чтобы заставить его "округлять до ближайшего целого", надо к числителю прибавить половинку знаменателя:
y0 = (j * 4096 + 67) / 135 - 16384;
x0 = -32768 + (i * 4096 + 67) / 135;
(кстати, деление - очень мерзопакостная операция ещё и в том, что округляет скорее не "вниз", а "к нулю", т.е в большинстве реализаций -8 / 9 также даст 0, что нарушает "изотропность" операций - вокруг нуля имеем "зону нечувствительности" расширенную вдвое. Здесь мы по счастью делим только положительные числа - одной проблемой меньше...)
Для осуществления на 32-разрядном процессоре, у которого есть команда деления - это уже неплохие формулы. Они дадут правильный результат, в пределах тех 16 бит, в которых мы храним результат.
Но когда у нас ПЛИС, то даже такие формулы кажутся слишком сложными: нам нужно делать умножение с сохранением результата в 32-битном аккумуляторе, затем ещё и провести деление - мягко говоря не самая любимая операция для реализации "в железе" (причём как ни странно, целочисленное деление оказывается в чём-то даже сложнее деления с плавающей запятой, поскольку там оба операнда фактически принимают значения от 1 до 2, а тут - какие угодно!).
По счастью, наша задача чуточку проще: мы знаем, что первоначально каждая из переменных i,j равна нулю, а затем мы прибавляем единичку к i во время движения по строке, после чего снова зануляем её и прибавляем единичку к j при переходе на новую строку.
Так что всё, что нам нужно - это правильно обновлять значения x0, y0 на каждом шаге.
Немножко перепишем выражения:
y0 = 30 * j + (46 * j + 67) / 135 - 16384;
x0 = 30 * i + (46 * j + 67) / 135 - 32768;
Мы видим, что на каждом шаге к переменной прибавляется 30 (за счёт первого слагаемого), и ещё должно добавляться 46 / 135 ≈ 0,34, но поскольку мы считаем в целых числах, то эти 0,34 будут "пропадать", но достаточно накопившись, дадут прибавление ещё единицы.
Собственно, если честно протабулировать значения y0, мы получаем:
-16384
-16354 (+30)
-16323 (+31)
-16293 (+30)
-16263 (+30)
-16232 (+31)
Получается, что на каждом шаге нам нужно лишь "принять решение" - прибавлять ли единицу? Само прибавление осуществить легко - в сумматоре есть вход переноса с младших разрядов - а ведь именно это у нас и происходит по сути :)
Чтобы принять решение, мы по сути храним "остаток от деления". Во время инициализации у нас уже осталась дробь 67/135, которую мы округлили до нуля. Знаменатель всегда будет 135 - его мы хранить не должны, он "сам собой подразумевается". А вот значение 67 мы запомним.
Пошёл первый шаг, и теперь дробь стала (67 + 46) / 135 = 113 / 135. Числитель всё ещё меньше знаменателя, поэтому единицу не прибавляем. Запоминаем значение 113.
На следующем шаге мы снова прибавляем к числителю 46, получая дробь (113 + 46) / 135 = 159 / 135. Числитель стал больше знаменателя - значит, мы прибавляем единицу к значению y0, а из нашего "остатка" эту единицу надо убрать, для чего мы вычитаем из числителя знаменатель, т.е получаем (159 - 135) / 135 = 24 / 135.
Вот так и работаем - ничего кроме сложений, вычитаний (а это одно и то же по сути) и сравнений нам не нужно. Есть целое семейство методов, схожих с этим, с помощью которых рисуют прямые, окружности, кривые Безье и прочие "графические примитивы" на растровом изображении.
Осталось только максимально "приспособить" эту задачу под специфику ПЛИС. Сравнение двух чисел "в лоб" - операция довольно неприятная, лучше было бы, если достаточно было посмотреть на знак числа. Для этого нужно "сдвинуть" дробь, чтобы именно "пересечение нуля" означало необходимость добавить единицу.
То есть, работа должна начинаться не с 67 / 135, а с 67 / 135 - 1 = -68 / 135. Прибавляем 46 первый раз - получаем -22 / 135 (пока ещё минус), прибавляем 46 второй раз - получаем 24 / 135 - это положительное число, поэтому прибавляем единицу к y0, вычитаем из числителя 135 и приходим к -111 / 135.
Вероятно, лучше всего будет ещё и "инвертировать" нашу дробь, чтобы старший разряд нашего "остатка" (в дополнительном коде появление единицы в старшем разряде означает - число стало отрицательным) напрямую подключался к входу переносу сумматора, не требуя дополнительной логики. То есть, начинаем мы с 68 / 135, на каждом шагу вычитаем 46, а если результат стал отрицательным - прибавляем единицу к y0, а к числителю прибавляем 135.
Наконец, замечаем, что не нужно спешить с прибавлением 135. Пусть по окончании шага результат остаётся отрицательным - просто в следующий раз, если числитель отрицательный, то вместо вычитания 46 мы будем прибавлять к нему 89. Это, к тому же, уменьшит диапазон чисел, который нам надо представлять в этом "остатке". Теперь самым большим числом может оказаться только 88 (ведь мы прибавляем значения лишь к отрицательным числам), а самым маленьким - число (-46), поскольку в худшем случае мы однажды придём к нулю, и в следующий раз вычтем из него 46. В нашем примере хватит 8-битного числа со знаком.
Пора приступать к кодированию на Verilog!
Первым делом, поиграемся с функциями и выражениями, исполняющимися на этапе компиляции:
module EqualPartsDivider (input clk, input sclr, input ce, output [15:0] Q);
parameter InitValue = -16384; //-16384
parameter Numerator = 32768;
parameter Denominator = 1080;
localparam CommonDivider = GCD (Numerator, Denominator);
localparam SimpleNum = Numerator / CommonDivider;
localparam SimpleDenom = Denominator / CommonDivider;
localparam WholePart = SimpleNum / SimpleDenom;
localparam Residue = SimpleNum - WholePart * SimpleDenom;
localparam RInitValue = (SimpleDenom + 1) / 2 - Residue;
localparam DenomMinusRes = SimpleDenom - Residue;
function integer GCD (input integer A, input integer B);
integer r;
begin
while (B != 0) begin
r = A % B;
A = B;
B = r;
end
GCD = A;
end
endfunction
initial begin
$display ("CommonDivider=%d", CommonDivider);
$display ("SimpleNum=%d", SimpleNum);
$display ("SimpleDenom=%d", SimpleDenom);
$display ("WholePart=%d", WholePart);
$display ("Residue=%d", Residue);
$display ("RInitValue=%d", RInitValue);
$display ("DenomMinusRes=%d", DenomMinusRes);
end
endmodule
Мы решили "довериться" компилятору в плане нахождения всех значений. Мы задаём ему 3 параметра: начальное значение, числитель и знаменатель, т.е говорим, что хотим получать значения выражения
Q = InitValue + Round (Numenator * i / Denominator);
для последовательных значений i.
Дальше это дело компилятора сообразить, где что надо прибавлять или вычитать.
Мы в кои-то веки написали функцию нахождения общего делителя (алгоритм Эвклида), и применять её мы будем только на этапе компиляции.
(вообще, функции могут быть доступны и в рантайме, но с этим связано множество проблем - иногда код просто не будет синтезироваться, поскольку не знает, сколько итераций будет внутри, а иногда даже синтезируется, но окажется чрезвычайно громоздким и медленным)
Затем мы сократили числитель и знаменатель, обозвав эти значения SimpleNum и SimpleDenom.
После этого мы выделили целую часть WholePart, и остаток Residue, и нашли исходное значение для нашего вспомогательного регистра, RInitValue (об этом чуточку позже). И напоследок, находим, что мы должны прибавить к остатку, когда он с прошлого раза остался отрицательным, это число мы назвали DenomMinusRes (знаменатель минус остаток).
Чтобы проверить правильность данных операций, мы решили задействовать несинтезируемую конструкцию $display - Квартус её тупо игнорирует, но есть симуляторы verilog, которые воспримут её и "напишут в консоль" то, что мы попросим. В том числе, есть онлайн-ресурсы, например, https://www.tutorialspoint.com/compile_verilog_online.php
(многие ютуберы показывают примеры на EDAboard, но там требуется регистрация, причём они ненавидят бесплатные почтовые сервисы и почтовые ящики, лежащие на серверах "конкурентов". Ну их нафиг...)
Вот что мы получаем:
$iverilog -o main *.v
$vvp main
CommonDivider= 8
SimpleNum= 4096
SimpleDenom= 135
WholePart= 30
Residue= 46
RInitValue= 22
DenomMinusRes= 89
135, 30, 46, 89 - эти числа мы совсем недавно использовали :)
Разберём заголовок модуля:
module EqualPartsDivider (input clk, input sclr, input ce, output reg [15:0] Q = 1'b0);
Интерфейс довольно простой: у нас есть вход тактовой частоты clk, вход для синхронного сброса sclr и вход разрешения счёта ce. Результирующее число мы получаем на выходе Q, причём оно устанавливается по положительному фронту clk во время подачи сигнала разрешения счёта ce, и сохраняет верное значение до последующих действий (сброса или счёта). Таким образом, "отлавливать", когда значение на выходе верное, нам не нужно.
Выход мы сразу определили как регистр. Ещё один регистр нужен для хранения "остатка". Его ширина может варьироваться в зависимости от дроби, которую мы пытаемся посчитать. Максимальное отрицательное значение не может превысить числа (-Residue), максимальное положительное значение - числа DenomMinusRes. Мы можем найти требуемую ширину регистра следующим образом:
localparam ResidueWidth1 = `CLOG2((DenomMinusRes + 1));
localparam ResidueWidth2 = `CLOG2(Residue);
localparam ResidueWidth = 1 + ((ResidueWidth1 > ResidueWidth2)? ResidueWidth1 : ResidueWidth2);
Макрос `CLOG2 мы сделали в части 6 - чрезвычайно полезная штука! Но буквально сейчас на ней обжёгся - поначалу у меня были строки
localparam ResidueWidth1 = `CLOG2(DenomMinusRes + 1) + 1;
localparam ResidueWidth2 = `CLOG2(Residue) + 1;
localparam ResidueWidth = (ResidueWidth1 > ResidueWidth2)? ResidueWidth1 : ResidueWidth2;
и "+1" в конце каждой строки у меня "сожралось", потому что `CLOG2 - это прямая подстановка выражения
(x <= 2) ? 1 : \
(x <= 4) ? 2 : \
(x <= 8) ? 3 : \
(x <= 16) ? 4 : \
(x <= 32) ? 5 : \
(x <= 64) ? 6 : \
(x <= 128)? 7 : \
(x <= 256)? 8 : \
(x <= 512)? 9 : \
(x <= 1024)? 10 : \
(x <= 2048)? 11 : \
(x <= 4096)? 12 : \
(x <= 8192)? 13 : \
(x <= 16384)? 14 : \
(x <= 32768)? 15 : \
(x <= 65536)? 16 : \
(x <= 'h20000)? 17 : \
(x <= 'h40000)? 18 : \
(x <= 'h80000)? 19 : \
(x <= 'h100000)? 20 : \
(x <= 'h200000)? 21 : \
(x <= 'h400000)? 22 : \
(x <= 'h800000)? 23 : \
(x <= 'h1000000)? 24 : \
(x <= 'h2000000)? 25 : \
(x <= 'h4000000)? 26 : \
(x <= 'h8000000)? 27 : \
(x <= 'h10000000)? 28 : \
(x <= 'h20000000)? 29 : \
(x <= 'h40000000)? 30 : \
(x <= 'h80000000)? 31 : \
(x <= 'h0FFFFFFFF)? 32 : \
-1
в строку. Единица прибавилась к последнему варианту. Чтобы такого не происходило, надо скобки ставить, чем больше - тем лучше! Всё как в препроцессоре C.
В нашем примере, ширина регистра составит 8 бит. Мы будем рассматривать его как число в дополнительном коде, т.е значения от -128 до +127.
Напишем, наконец, логику нашего модуля! На "чистом верилоге" она очень проста:
reg [ResidueWidth - 1 : 0] D = 1'b0; //Decision
wire IsNeg = D[ResidueWidth - 1];
assign DD = D;
always @(posedge clk) if (ce | sclr) begin
D <= sclr? RInitValue : IsNeg? D + DenomMinusRes : D - Residue;
Q <= sclr? InitValue : Q + WholePart + IsNeg;
end
И действительно, этот код заработает так, как нам надо. Но результат синтеза меня не радует: для наших параметров получилось 78 ЛЭ. Точнее, само по себе число 78 - не такое уж страшное (если только это не заголовок Zlib), хуже другое: в логах мы видим, что в этих двух строках компилятор усмотрел аж 4 сумматора!
Как всегда, он не понял, что прибавление единички не требует дополнительного сумматора, ведь у нас есть вход переноса, ну ё-моё! И не понял, что можно применить один сумматор в первой строке, просто подставить ему разные значения.
Так что нам ничего не остаётся, как вставить сумматоры ручками:
reg [ResidueWidth - 1 : 0] D = 1'b0; //Decision
assign DD = D;
wire IsNeg = D[ResidueWidth - 1];
wire [ResidueWidth - 1 : 0] ValueToAdd = IsNeg? DenomMinusRes : (-Residue);
wire [ResidueWidth - 1: 0 ] DResult;
lpm_add_sub ResidueAdder ( .dataa (D),
.datab (ValueToAdd),
.result (DResult) );
defparam
ResidueAdder.lpm_direction = "ADD",
ResidueAdder.lpm_hint = "ONE_INPUT_IS_CONSTANT=NO,CIN_USED=NO",
ResidueAdder.lpm_representation = "SIGNED",
ResidueAdder.lpm_type = "LPM_ADD_SUB",
ResidueAdder.lpm_width = ResidueWidth;
wire [15:0] QResult;
lpm_add_sub MainAdder ( .dataa (Q),
.datab (WholePart),
.cin (IsNeg),
.result (QResult) );
defparam
MainAdder.lpm_direction = "ADD",
MainAdder.lpm_hint = "ONE_INPUT_IS_CONSTANT=YES,CIN_USED=YES",
MainAdder.lpm_representation = "SIGNED",
MainAdder.lpm_type = "LPM_ADD_SUB",
MainAdder.lpm_width = 16;
always @(posedge clk) if (ce | sclr) begin
D <= sclr? RInitValue : DResult;
Q <= sclr? InitValue : QResult;
end
Теперь мы получаем те же самые результаты, но используя лишь 50 ЛЭ, а что даже важнее - мы вдвое сократили комбинаторные пути, так что этот модуль спокойненько компилится и на медленном кристалле "-3", при заданной частоте в 80 МГц (даже остался запас по частоте аж в 10 МГц). Исходная схема уже ругалась...
Под спойлером - код модуля целиком (включая заголовок и толпу локальных параметров):
[Spoiler (click to open)]
`include "math.v"
module EqualPartsDivider (input clk, input sclr, input ce, output reg [15:0] Q = 1'b0, output [7:0] DD);
parameter InitValue = -16384; //-16384
parameter Numerator = 32768;
parameter Denominator = 1080;
localparam CommonDivider = GCD (Numerator, Denominator);
localparam SimpleNum = Numerator / CommonDivider;
localparam SimpleDenom = Denominator / CommonDivider;
localparam WholePart = SimpleNum / SimpleDenom;
localparam Residue = SimpleNum - WholePart * SimpleDenom;
localparam RInitValue = (SimpleDenom + 1) / 2 - Residue;
localparam DenomMinusRes = SimpleDenom - Residue;
localparam ResidueWidth1 = `CLOG2(DenomMinusRes + 1);
localparam ResidueWidth2 = `CLOG2(Residue);
localparam ResidueWidth = 1 + ((ResidueWidth1 > ResidueWidth2)? ResidueWidth1 : ResidueWidth2);
reg [ResidueWidth - 1 : 0] D = 1'b0; //Decision
assign DD = D;
wire IsNeg = D[ResidueWidth - 1];
wire [ResidueWidth - 1 : 0] ValueToAdd = IsNeg? DenomMinusRes : (-Residue);
wire [ResidueWidth - 1: 0 ] DResult;
lpm_add_sub ResidueAdder ( .dataa (D),
.datab (ValueToAdd),
.result (DResult) );
defparam
ResidueAdder.lpm_direction = "ADD",
ResidueAdder.lpm_hint = "ONE_INPUT_IS_CONSTANT=NO,CIN_USED=NO",
ResidueAdder.lpm_representation = "SIGNED",
ResidueAdder.lpm_type = "LPM_ADD_SUB",
ResidueAdder.lpm_width = ResidueWidth;
wire [15:0] QResult;
lpm_add_sub MainAdder ( .dataa (Q),
.datab (WholePart),
.cin (IsNeg),
.result (QResult) );
defparam
MainAdder.lpm_direction = "ADD",
MainAdder.lpm_hint = "ONE_INPUT_IS_CONSTANT=YES,CIN_USED=YES",
MainAdder.lpm_representation = "SIGNED",
MainAdder.lpm_type = "LPM_ADD_SUB",
MainAdder.lpm_width = 16;
always @(posedge clk) if (ce | sclr) begin
D <= sclr? RInitValue : DResult;
Q <= sclr? InitValue : QResult;
end
function integer GCD (input integer A, input integer B);
integer r;
begin
while (B != 0) begin
r = A % B;
A = B;
B = r;
end
GCD = A;
end
endfunction
initial begin
$display ("CommonDivider=%d", CommonDivider);
$display ("SimpleNum=%d", SimpleNum);
$display ("SimpleDenom=%d", SimpleDenom);
$display ("WholePart=%d", WholePart);
$display ("Residue=%d", Residue);
$display ("RInitValue=%d", RInitValue);
$display ("DenomMinusRes=%d", DenomMinusRes);
$display ("ResidueWidth=%d", ResidueWidth);
end
endmodule
Ну и объясним напоследок значение RInitValue - оно сюда внесено "с опережением на шаг"! Т.е в наших размышлениях мы начинали со значения 68 (половинка от знаменателя, округлённая "вверх"), а здесь сразу же вычли ещё 46, получив значение 22. Это позволяет сумматору для Q брать значение не с выхода сумматора для D (как мы делали в своих устных вычислениях - сначала вычитаем здесь, потом, если получили отрицательное значение - прибавим здесь единицу), а с выхода регистра. А для этого нужно, чтобы к началу такта в регистре уже хранилось правильное значение. По счастью, в этой задаче мы заранее знали, как его посчитать - вот и посчитали заблаговременно!
Для нашего генератора множества Мандельброта нужно два таких модуля, с разными параметрами. Заметим, что они не только вычисляют выражение, но и служат регистрами для хранения текущего Z0.
Смотрите в следующей части - "аппаратная" реализация выражения Z = Z2 + Z0, с проверкой на переполнение.
Poll Целочисленная арифметика