Домашка по арифметике

в 22:46, , рубрики: fpga, Алгоритмы, математика, реализация умножения

Лёшенька, Лёшенька, сделай одолжение!
Выучи, Алёшенька, таблицу умножения !

Агния Барто

Домашка по арифметике - 1

Сначала задачка для первоклассника. Дано некоторое положительное число. Нужно умножить на него другое число, заранее неизвестное. Вопрос, как посоветуют это сделать благородные доны ??? Бывалый разраб наверняка скажет, мол мужик, ставь умножитель и не парь мне мОзги. И возможно будет в корне неправ! Ибо кроме монстров от Alterra и Xilinx существует ещё и такое замечательное семейство как iCE-40 от Lattice. Ультрамикропотребляющее. Очень дешевое. Да вот беда, больно мелкие они, и увы, умножителей там нет. Я столкнулся с этим года 4 назад, когда портировал некий ADPCM-кодек с ассемблера adsp-2185 на такой кристалл.

Нет, без умножителя общего назначения(для двух переменных) там разумеется не обошлось. Пришлось сделать его ручками и он занял половину кристалла. Беда в том, что молотил он у меня в каждом такте, т.е. вообще без свободных окон. А кроме умножений переменных в алгоритме были фиксированные коэффициенты, на которые тоже надо было умножать. И использовать для этого умножитель не было никакой возможности, просто не оставалось свободного окна. А поскольку борьба в проекте шла за каждую ячейку, разумеется меня очень интересовало, как бы извернуться, и используя свойства этих коэффициентов сделать оптимальную схему умножения на них(не умножитель общего назначения, а устройство умножения на заданное число !). Ну и как полная сбыча мечт — чтобы схемы умножения на разные коэффициенты максимально использовали какие-то общие ресурсы кристалла.

Итак, от задачки для первоклассника мы плавно перешли к задаче для инженера и математика. Как для заданного конкретного числа построить оптимальную схему умножения на него?

Вспомним немного, как мы вообще умножаем число А на В.

  1. Для начала положим результат равным нулю.
  2. Пробежимся по разрядам A.
  3. Если в разряде 0, ничего не делаем.
  4. Если в разряде 1, сдвигаем B на соответствующее число разрядов и прибавляем к результату.

Итого, чем меньше в разрядах А единиц, тем проще на него умножать. Умножить на 2, 4 или 8 совсем легко. Нужен только сдвиг. Немного труднее умножить на 3, 5 или 10. Потребуются одно сложение. Ещё труднее умножить на 11. Чем больше единиц, тем труднее.

Ну а как насчёт умножения на 255? Целых 8 единиц, кошмарная задача, правда? А вот и фигушки! Сделаем хитрый финт ушами. Умножим наше В на 256 (т.е. сдвинем на 8 разрядов) и вычтем его же из результата. Получается несмотря на грозный вид, умножить на 255 столь же легко, как и на 10. Столь же легко умножать и на 254, и на 240, и на 224 (если хотите, проверьте !).

Итого, вычитание помогает умножению. Запомним эту ценную мысль. А пока назовём трудностью числа минимальное число операций сложения (или вычитания), необходимых для умножения на него. Сдвиги не учитываем. Ибо в контексте задачи(у нас схемотехника, а не программирование !) они получаются даром. Сформулируем пока одно очевидное, но полезное утверждение:

Трудность числа меньше или равна количеству единиц в его двоичной записи.
В самом деле, если мы без фокусов с вычитанием честно умножим, как в учили в школе, число сложений будет равно числу единиц в двоичной записи A. А если начнём хитрить и к тому же нам повезёт, то сократим число операций.

Интересное следствие из этого утверждения:
Трудность любого n-разрядного числа всегда строго меньше n.
Ведь наибольшее число единиц в n-разрядном числе равно n, у числа подобного 255. С другой стороны понятно, что фокус, как с умножением на 255, мы можем повторить для чисел любой разрядности. Это открывает перспективы и для оптимизации умножителей общего назначения.

Придадим нашим рассуждениям более регулярный вид. Для начала пусть у нас уже есть какая-то схема умножения на А. Возьмём в качестве В единицу. Тогда все слагаемые сдвинутся, сложатся и вычтутся таким образом, что в результате получится то же самое А. Итого, наша схема умножения на A есть ни что иное, как представление числа A в виде суммы степеней двойки, в котором слагаемые могут быть и со знаком + и со знаком -. А можно сгруппировать вместе положительные и отрицательные слагаемые, и сказать что схема описывается разностью неких чисел A+ и A-, равной A. Это плохо… Члены разности могут быть как угодно велики. Подумаем, есть ли какие-то соображения, позволяющие их ограничить?
Заметим прежде всего, что каждая степень двойки может входить в схему только один раз. В самом деле, если она входит дважды с одним знаком, это можно заменить степенью двойки на единицу бОльшей. Если она входит дважды с разными знаками, они взаимно вычитаются. Итого, схема очень напоминает обычное двоичное представление A. С той лишь разницей, что его «биты»(в дальнейшем будем называть их триты), могут принимать не два, а три значения: 1, 0 и -1.

Итого. Если A есть n-разрядное двоичное число, схема умножения на него описывается суммой:

$A={t_m}2^m + {t_{m-1}}2^{m-1} + ....+ {t_1}2 + {t_0}$

где ${t_i}$ — триты, принимающие значения 1, 0 и -1, а m некое пока неизвестное число.

В дальнейшем будем называть такое выражение троичной схемой или просто схемой или схемой умножения для числа A.

Попытаемся найти m. Как видно из примера с умножением на 255, m никак не меньше n. Посмотрим, может ли оно быть равно n+1.
Пускай A как и прежде, положительное n-разрядное число. Тогда схема умножения задаётся как:

$A={t_{n+1}}2^{n+1} + {t_{n}}2^{n} + {t_{n-1}}2^{n-1} ....+ {t_1}2 + {t_0}$

Вспомним, что $2^k + 2^{k-1} +....+ 2 + 1=2^{k+1}-1$
Поэтому старший трит никак не может быть равен -1. Ведь даже если все остальные будут 1, сумма всё равно будет равна -1. А по условию A положительно. Пусть теперь старший трит равен 1. Но в этом случае следующий трит должен быть равен -1. Иначе сумма получится слишком большой. В самом деле, если следующий трит равен 0, сумма получится не меньше чем $2^{n+1} - (2^{n-1}+...+1)=2^{n+1} - (2^n - 1)=2^n+1$
В то время как n-разрядное A не больше $2^n - 1$. Но если старший трит равен 1, а следующий -1, то сумма двух старших членов равна $2^{n+1}-2^n=2^n$. А значит старший трит равен 0.

Таким образом доказано важное утверждение: m=n.

Иными словами для представления любой схемы умножения на n-разрядное A достаточно n+1 степеней двойки — от 0 до n(на одну больше, чем для двоичного представления), что сразу избавляет нас от бесконечности.

Теперь можно попытаться вычислить трудности(см выше) каких-то конкретных чисел. Самое простое это поступить так, как мы поступаем выписывая по порядку двоичные числа. Только во-первых у нас тут не биты а триты (три возможных значения), во-вторых некоторые комбинации тритов дают отрицательные числа и их надо отбрасывать, в третьих некоторые комбинации могут давать совпадающие числа. Это означает что для такого числа существует несколько схем.

Писать будем на питоне. Не потому что я такой его фанат, а потому что с ним крайне удобно работать в блокнотах jupyter. Для начала напишем класс TriScheme, работающий с введёнными выше троичными схемами, функцию allSchemes, вычисляющую с его помощью все схемы для чисел заданной разрядности простым перебором, и функцию printSchemes, распечатывающую результат:

TriScheme

class TriScheme():     
    def __init__(self, buf:list): # создает из массива
        self.buf=[]
        for i in range(0, len(buf)):
            s=0
            if (type(buf[i]) == int) or (type(buf[i]) == float):
                if buf[i]   > 0: s=1
                elif buf[i] < 0: s=-1
            self.buf.append(s)
            
    def zero(self): # обнуляет
        self.buf=[0]*len(self.buf)
        
    def __repr__(self): # представление в виде строки
        s=""
        for i in range(1, len(self.buf)+1):            
            if   self.buf[-i]==1:  s=s+"+"
            elif self.buf[-i]==0:  s=s+"0"
            elif self.buf[-i]==-1: s=s+"-"
        return s
            
    def inc(self):  # инкремент (переход к следующей)
        for i in range(0, len(self.buf)):
            if (self.buf[i]+1) < 2:
                self.buf[i]+=1
                break
            for j in range(0, i+1):
                self.buf[i]=-1
                
    def __int__(self):  # выдаёт значение
        m=1
        s=0
        for i in range(0, len(self.buf)):
            s+=self.buf[i]*m
            m=m*2
        return s
    
    def isMax(self): # проверяет что схема максимальна
        s=0
        for i in range(0, len(self.buf)):
            s+=self.buf[i]
        return (s == len(self.buf))
    
    def isMaxDig(self): # проверяет, что схема максимальная, соответствующая n-разрядному числу
        s=[0]*len(self.buf)
        s[0]=-1
        s[-1]=1
        return (self.buf == s)
        
    def ones(self): # выдает число ненулевых элементов (трудность схемы)
        s=0
        for i in range(0, len(self.buf)):
            if self.buf[i] != 0:
                s+=1
        return s
    
    def minPos(self):  # выдаёт номер минимальной ненулевой позиции
        for i in range(0, len(self.buf)):
            if self.buf[i] != 0:
                return i
                         
def allSchemes(dig): # считает все схемы для всех чисел разрядности dig
    sch=[]
    for i in range(0, 2**dig): sch.append([])
    ts=TriScheme([0]*(dig+1))    
    while True:
        sch[int(ts)].append(TriScheme(ts.buf))
        if ts.isMaxDig():
            break
        ts.inc()
    return sch

def printSchemes(sch):  # печатает список схем в формате число трудность количество_схем список_схем
    for i in range(0, len(sch)):
        s=sch[i]
        a=[]
        for j in range(0, len(s)): 
            a.append(s[j].ones())
        m=min(a)
        print(i, m, len(s), s)

# Вычисляем схемы для всех 4-разрядных чисел        
sch4=allSchemes(4)
printSchemes(sch4)

Вычислим схемы для всех 4-разрядных чисел:

0 0 1 [00000]
1 1 5 [0000+, 000+-, 00+--, 0+---, +----]
2 1 4 [000+0, 00+-0, 0+--0, +---0]
3 2 7 [000++, 00+-+, 00+0-, 0+--+, 0+-0-, +---+, +--0-]
4 1 3 [00+00, 0+-00, +--00]
5 2 8 [00+0+, 00++-, 0+-0+, 0+-+-, 0+0--, +--0+, +--+-, +-0--]
6 2 5 [00++0, 0+-+0, 0+0-0, +--+0, +-0-0]
7 2 7 [00+++, 0+-++, 0+0-+, 0+00-, +--++, +-0-+, +-00-]
8 1 2 [0+000, +-000]
9 2 7 [0+00+, 0+0+-, 0++--, +-00+, +-0+-, +-+--, +0---]
10 2 5 [0+0+0, 0++-0, +-0+0, +-+-0, +0--0]
11 3 8 [0+0++, 0++-+, 0++0-, +-0++, +-+-+, +-+0-, +0--+, +0-0-]
12 2 3 [0++00, +-+00, +0-00]
13 3 7 [0++0+, 0+++-, +-+0+, +-++-, +0-0+, +0-+-, +00--]
14 2 4 [0+++0, +-++0, +0-+0, +00-0]
15 2 5 [0++++, +-+++, +0-++, +00-+, +000-]

Здесь в начале строки число, за ним его трудность (минимальное количество ненулевых элементов в схеме), затем число схем, и наконец список схем.

В схеме + обозначает 1, 0 — 0, и - — -1. Старший трит слева (как при обычном написании чисел).

Из таблицы видно во-первых, что только числа 13 и 11 имеют трудность 3 (максимальную для 4-разрядных чисел, как доказано выше). И во-вторых, что число различных схем для каждого числа довольно велико. Это говорит во-первых о том, что умножение чаще всего операция сильно более простая чем нас учили в школе. И во-вторых о том, что для реализации устройств умножения на несколько констант с использованием общих ресурсов кристалла, выбор схем достаточно велик. Ещё интереснее выглядят данные по более длинным числам. Так для 14-разрядных чисел максимальная трудность получается равной 8. А максимальное число схем — 937.

Всё бы хорошо, но приведённый выше алгоритм слишком затратен. Его сложность растёт как $3^n$ в зависимости от разрядности чисел n. Для 8 разрядов считается мгновенно. Для 14 — минуты. Для 16 — часы. Если нужно считать схемы ДЛЯ ВСЕХ n-разрядных чисел, увы, больше ничего сделать нельзя, кроме как переписать его на С и взять компьютер помощнее. Впрочем алгоритм прекрасно распараллеливается и наверняка может быть запущен на GPU или на кластере.

Для проектирования конкретных железок как правило требуются списки схем для конкретных чисел. Причём желательно ВСЕ, а не только оптимальные. Ибо какую выбрать, может зависеть от обстоятельств(например устройство умножения на несколько констант). И вот тут можно предложить алгоритм гораздо более гуманный. Ещё раз вспомним замечательное равенство: $2^k + 2^{k-1} +....+ 2 + 1=2^{k+1}-1$.

В контексте задачи это означает следующее. Пусть известны несколько старших тритов схемы. Все младшие триты, начиная с позиции k неизвестны и предварительно считаются равными нулю. Тогда изменяя каким угодно образом неизвестные триты, мы изменим значение схемы не более чем на плюс/минус $2^{k+1}-1$. Причем с продвижением к младшим тритам величина этой неопределённости уменьшается. Это позволяет искать схемы методом последовательных приближений, трит за тритом.

Пусть дано положительное число A. Надо найти все схемы для него среди n-разрядных чисел. Абсолютную величину разности между A и значением некоторой схемы, назовём погрешностью схемы. Возьмем для начала n+1-разрядную схему, описывающую n-разрядные числа, все триты которой неизвестны, и изначально положены равными нулю. И начнём определять триты, один за другим, начиная со старшего. В k-й позиции схема может породить три новых схемы — со значениями k-го трита 1, 0 и -1. Оставим в списке схем для продолжения только те из них, погрешность которых не превышает $2^k-1$. С получившимся списком схем перейдём к шагу в позиции k-1. Таким образом, в результирующем списке (в позиции 0) будут ВСЕ возможные схемы, значение которых равно A. Давайте напишем такую функцию calcSchemes.

calcSchemes

def calcSchemes(a, n=-1):
    m=1
    nn=1
    while m < a:
        nn+=1
        m=m*2
    if n < nn:
        n=nn
    sch=[TriScheme([0]*n)]
    for i in range(0, n):
        tmp=[]
        pos=n-i-1
        for j in range(0, len(sch)):
            for k in range(-1, 2):
                ts=sch[j]
                ts.buf[pos]=k
                d=abs(a - int(ts))
                if d <= (2**pos - 1):
                    tmp.append(TriScheme(ts.buf))
        sch=tmp   
    return sch 

Ну и наконец обещанная сбыча мечт.

В том проекте было два коэффициента, на которые требовалось умножать: 16351 и 16318. Используя calcSchemes, найдём схемы для этих коэффициентов:

Схемы

Для 16351
0++++++++0+++++
0+++++++++-++++
0+++++++++0-+++
0+++++++++00-++
0+++++++++000-+
0+++++++++0000-
+-+++++++0+++++
+-++++++++-++++
+-++++++++0-+++
+-++++++++00-++
+-++++++++000-+
+-++++++++0000-
+0-++++++0+++++
+0-+++++++-++++
+0-+++++++0-+++
+0-+++++++00-++
+0-+++++++000-+
+0-+++++++0000-
+00-+++++0+++++
+00-++++++-++++
+00-++++++0-+++
+00-++++++00-++
+00-++++++000-+
+00-++++++0000-
+000-++++0+++++
+000-+++++-++++
+000-+++++0-+++
+000-+++++00-++
+000-+++++000-+
+000-+++++0000-
+0000-+++0+++++
+0000-++++-++++
+0000-++++0-+++
+0000-++++00-++
+0000-++++000-+
+0000-++++0000-
+00000-++0+++++
+00000-+++-++++
+00000-+++0-+++
+00000-+++00-++
+00000-+++000-+
+00000-+++0000-
+000000-+0+++++
+000000-++-++++
+000000-++0-+++
+000000-++00-++
+000000-++000-+
+000000-++0000-
+0000000-0+++++
+0000000-+-++++
+0000000-+0-+++
+0000000-+00-++
+0000000-+000-+
+0000000-+0000-
+00000000--++++
+00000000-0-+++
+00000000-00-++
+00000000-000-+
+00000000-0000-

Для 16318
0+++++++0+++++0
0++++++++-++++0
0++++++++0-+++0
0++++++++00-++0
0++++++++000-+0
0++++++++0000-0
+-++++++0+++++0
+-+++++++-++++0
+-+++++++0-+++0
+-+++++++00-++0
+-+++++++000-+0
+-+++++++0000-0
+0-+++++0+++++0
+0-++++++-++++0
+0-++++++0-+++0
+0-++++++00-++0
+0-++++++000-+0
+0-++++++0000-0
+00-++++0+++++0
+00-+++++-++++0
+00-+++++0-+++0
+00-+++++00-++0
+00-+++++000-+0
+00-+++++0000-0
+000-+++0+++++0
+000-++++-++++0
+000-++++0-+++0
+000-++++00-++0
+000-++++000-+0
+000-++++0000-0
+0000-++0+++++0
+0000-+++-++++0
+0000-+++0-+++0
+0000-+++00-++0
+0000-+++000-+0
+0000-+++0000-0
+00000-+0+++++0
+00000-++-++++0
+00000-++0-+++0
+00000-++00-++0
+00000-++000-+0
+00000-++0000-0
+000000-0+++++0
+000000-+-++++0
+000000-+0-+++0
+000000-+00-++0
+000000-+000-+0
+000000-+0000-0
+0000000--++++0
+0000000-0-+++0
+0000000-00-++0
+0000000-000-+0
+0000000-0000-0

Нам повезло! Оба коэффициента имеют трудность 3. Выбираем для обоих оптимальные схемы:
Для 16318 +0000000-0000-0 и для 16351 — +00000000-0000-. Нам ещё раз повезло! Обратите внимание на хвосты схем. Они для 16318 и 16351 отличаются только сдвигом влево на одну позицию. А значит устройство умножения на 16318 и 16351, включает только один добавочный мультиплексор, переключающий сдвинутый и не сдвинутый операнд на входе сумматора.

Давайте посмотрим результат. Напишем на верилоге три варианта устройства умножения на 16318 и 16351:

  1. «Школьный» вариант (только сложения и сдвиги)
  2. Оптимальная схема
  3. Оптимальная схема, использующая общие ресурсы

Для всех трёх вариантов выполним синтез, и посмотрим сколько на каждый из них будет затрачено ресурсов кристалла.

Verilog

/*-----------------------Школьный вариант----------------------*/
module mul_163xx(
        input[15:0] di,
	input s,
	output[31:0] dq
           );
reg[31:0] dd;
always @*
	if(s) dd= (di<<13)+(di<<12)+(di<<11)+(di<<10)+(di<<9)+(di<<8)+(di<<7)+(di<<4)+(di<<3)+(di<<2)+(di<<1) + (di<<6)+di;
	else  dd=(di<<13)+(di<<12)+(di<<11)+(di<<10)+(di<<9)+(di<<8)+(di<<7)+(di<<4)+(di<<3)+(di<<2)+(di<<1) + (di<<5); 
assign dq=dd;
endmodule

/*---------------------------Оптимальная схема -------------------------*/
module mul_163xx(
        input[15:0] di,
	input s,
	output[31:0] dq
           );
reg[31:0] dd;
always @*
	if(s) dd= (di<<14) - (di<<5) - di;
	else dd=(di<<14) - (di<<6) - (di<<1);
assign dq=dd; 
endmodule

/*--------------Оптимальная схема с общими ресурсами--------------*/
module mul_163xx(
        input[15:0] di,
	input s,
	output[31:0] dq
           );
wire[31:0] tail = (di<<5) + di;
reg[31:0] dd;
always @*
	if(s) dd= (di<<14) - tail;
	else dd=(di<<14) - (tail<<1);
assign dq=dd; 
endmodule

/*-------------------------------Тестбенч------------------------------*/
module mult_tb();

reg clk = 0;
always #100 clk = ~clk;

reg[15:0] di;
wire[31:0] dq;
reg s; 

mul_163xx mul (
	.di (di),
	.dq (dq),
	.s  (s)
         );

initial 
	begin
	clk=0;

	s=0;
	$display("s=0");
	
	di=1;
	@(posedge clk);
	$display("%d", dq);

	di=10;
	@(posedge clk);
	$display("%d", dq);

	di=100;
	@(posedge clk);
	$display("%d", dq);

	di=1000;
	@(posedge clk);
	$display("%d", dq);

	s=1;
	$display("s=1");

	di=1;
	@(posedge clk);
	$display("%d", dq);

	di=10;
	@(posedge clk);
	$display("%d", dq);

	di=100;
	@(posedge clk);
	$display("%d", dq);

	di=1000;
	@(posedge clk);
	$display("%d", dq);

  	$finish;
	end

endmodule

/*------------------------------PCF-файл------------------------------*/
set_io  dq[0]  A1
set_io  dq[1]  A2
set_io  dq[2]  P16
set_io  dq[3]  M13
set_io  dq[4]  A5
set_io  dq[5]  A6
set_io  dq[6]  A7
set_io  dq[7]  L13
set_io  dq[8]  A9
set_io  dq[9]  A10
set_io  dq[10]  A11
set_io  dq[11]  M14
set_io  dq[12]  P15
set_io  dq[13]  N16
set_io  dq[14]  A15
set_io  dq[15]  A16
set_io  dq[16]  B1
set_io  dq[17]  B2
set_io  dq[18]  B3
set_io  dq[19]  B4
set_io  dq[20]  B5
set_io  dq[21]  B6
set_io  dq[22]  B7
set_io  dq[23]  B8
set_io  dq[24]  B9
set_io  dq[25]  B10
set_io  dq[26]  B11
set_io  dq[27]  B12
set_io  dq[28]  B13
set_io  dq[29]  B14
set_io  dq[30]  B15
set_io  dq[31]  B16
set_io  di[0]  C1
set_io  di[1]  C2
set_io  di[2]  C3
set_io  di[3]  C4
set_io  di[4]  C5
set_io  di[5]  C6
set_io  di[6]  C7
set_io  di[7]  C8
set_io  di[8]  C9
set_io  di[9]  C10
set_io  di[10]  C11
set_io  di[11]  C12
set_io  di[12]  C13
set_io  di[13]  C14
set_io  di[14]  D4
set_io  di[15]  C16
set_io  s  D1

Все три варианта работают правильно, умножая при 0 на входе s на 16318, а при 1 — на 16351. При этом yosys даёт для школьного варианта 488 ячеек, для оптимального — 206, и для варианта с общими ресурсами 202. Наверно он сам тоже что-то оптимизирует, хотя разница в 4 ячейки всё-таки есть. Как видим, разница со школьным вариантом очень приличная.

Ну и наконец. Возможно кому-то покажется излишним городить такой огород, только для того, чтобы элементарно сообразить, что 16318=16384-64-2, а 16351=16384-32-1. Однако во-первых числа могут быть и посложнее. Во-вторых, если устройство должно умножать несколько чисел, далеко не очевидно, что следует брать оптимальные схемы. Мне в том проекте просто повезло. В общем случае программа поиска схем может сильно помочь. Надеюсь статья кому-то оказалось полезной. И тот кто её прочитал, надеюсь не будет паниковать, если надо умножать, а умножителя нет.

Автор: Калягин Евгений Игоревич

Источник


* - обязательные к заполнению поля


https://ajax.googleapis.com/ajax/libs/jquery/3.4.1/jquery.min.js