БПФ - перевод с С++ на Делфи

[Семионн]

Буду очень благодарен тому, кто поможет исправить ошибки в моей попытке организовать умножение длинных чисел с Быстрым Преобразованием Фурье. Само преобразование понимаю лишь отчасти, хотя перечитал литературы... От Википедии до Кнута) В итоге пробовал перевести этот алгоритм на С++ на Делфи. Подумал, что среди людей знающих плюсы больше тех кто знает паскаль, чем наоборот) Поэтому тема здесь.
Заранее благодарен!!

int rev[MAXN];
base wlen_pw[MAXN];
 
void fft (base a[], int n, bool invert) {
	for (int i=0; i<n; ++i)
		if (i < rev[i])
			swap (a[i], a[rev[i]]);
 
	for (int len=2; len<=n; len<<=1) {
		double ang = 2*PI/len * (invert?-1:+1);
		int len2 = len>>1;
 
		base wlen (cos(ang), sin(ang));
		wlen_pw[0] = base (1, 0);
		for (int i=1; i<len2; ++i)
			wlen_pw[i] = wlen_pw[i-1] * wlen;
 
		for (int i=0; i<n; i+=len) {
			base t,
				*pu = a+i,
				*pv = a+i+len2, 
				*pu_end = a+i+len2,
				*pw = wlen_pw;
			for (; pu!=pu_end; ++pu, ++pv, ++pw) {
				t = *pv * *pw;
				*pv = *pu - t;
				*pu += t;
			}
		}
	}
 
	if (invert)
		for (int i=0; i<n; ++i)
			a[i] /= n;
}
 
void calc_rev (int n, int log_n) {
	for (int i=0; i<n; ++i) {
		rev[i] = 0;
		for (int j=0; j<log_n; ++j)
			if (i & (1<<j))
				rev[i] |= 1<<(log_n-1-j);
	}
}

void multiply (const vector<int> & a, const vector<int> & b, vector<int> & res) {
	vector<base> fa (a.begin(), a.end()),  fb (b.begin(), b.end());
	size_t n = 1;
	while (n < max (a.size(), b.size()))  n <<= 1;
	n <<= 1;
	fa.resize (n),  fb.resize (n);
 
	fft (fa, false),  fft (fb, false);
	for (size_t i=0; i<n; ++i)
		fa[i] *= fb[i];
	fft (fa, true);
 
	res.resize (n);
	for (size_t i=0; i<n; ++i)
		res[i] = int (fa[i].real() + 0.5);int carry = 0;
	for (size_t i=0; i<n; ++i) {
		res[i] += carry;
		carry = res[i] / 10;
		res[i] %= 10;
	}
}

(брал с сайта: http://e-maxx.ru/algo/fft_multiply)
Мой выводящий нули вариант на паскале... а иногда вообще вылетающий в рантайм(

type complex = record
     re,im:double;
  end;
  type arcomp = array of complex;
  type arint = array of integer;

var a,b:arcomp;
r:arint;
s:string;
i,n1,n2,t:integer;
rev: arint;
wlen_pw:arcomp;

function umnoj(a,b:complex):complex;
begin
  umnoj.re:=a.re*b.re-a.im*b.im;
  umnoj.im:=a.re*b.im+a.im*b.re;
end;

function summ(a,b:complex):complex;
begin
  summ.re:=a.re+b.re;
  summ.im:=a.im+b.im;
end;

function diff(a,b:complex):complex;
begin
  diff.re:=a.re-b.re;
  diff.im:=a.im-b.im;
end;

function base(a,b:double):complex;
begin
  base.re:=a;
  base.im:=b;
end;

procedure fft (a:arcomp;invert:boolean; n:integer );
var len,len2,k,j,d:integer;
ang:double;
wlen,t:complex;
wlen_pw:arcomp;
procedure swap(var a,b:complex);
var c:complex;
begin
  c:=a;
  a:=b;
  b:=c;
end;
begin
	for i:=0 to n do
		if i<rev[i] then
			swap (a[i], a[rev[i]]);
  len:=2;
  setlength(wlen_pw,n);
	while (len<=n)do begin
    if invert then
		  ang := (-2)*PI/len
    else
		  ang := 2*PI/len;
		len2 := len shr 1;
		wlen.re:=cos(ang);
    wlen.im:=sin(ang);
		wlen_pw[0]:= base(1, 0);

		for i:=1 to len2-1 do
			wlen_pw[i] := umnoj(wlen_pw[i-1],wlen);
    i:=0;
    while i<n do begin
      k:=0; j:=0; d:=0;
      while k<>len2 do begin
        t:=umnoj(a[i+len2+j],wlen_pw[d]);
        a[i+len2+j]:=diff(a[i+k],t);
        inc(k); inc(j); inc(d);
      end;
      i:=i+len;
		end;
    len:=len shl 1;
	end;
	if invert then
		for i:=0 to n do begin
			a[i].re := a[i].re/n;
      a[i].im := a[i].im/n;
    end;
end;


procedure calc_rev (n,log_n:integer);
var i,j:integer;
begin
	for i:=0 to n-1 do begin
		rev[i]:=0;
		for j:=0 to log_n-1 do
			if (i and (1 shl j))<>0 then
				rev[i]:=rev[i]or (1 shl (log_n-1-j));
	end;
end;

procedure multiply (a,b:arcomp; var res: arint;len1,len2:integer);
var n,_as,bs:integer;
i,carry,k:integer;
fa,fb:arcomp;
function max(a,b:integer):integer;
begin
  if a>b then max:=a else max:=b;
end;

begin
  n:=1;
  setlength(fa,len1);
  setlength(fb,len2);
  for i:=0 to len1-1 do
    fa[i]:=a[i];
  for i:=0 to len2-1 do
    fb[i]:=b[i];
  _as:=len1;
  bs:=len2;
  while (n<max(_as,bs)) do
    n:=n*2;
  n:=n*2;
  setlength(fa,n);
  setlength(fb,n);
  setlength(rev,n);
  calc_rev(n,round(ln(n)));
  fft(fb,false,n);
  fft(fa,false,n);
  for i:= 0 to n-1 do
    fa[i]:= umnoj(fa[i],fb[i]);
  fft(fa,true,n);
  setlength(res,n);
  for i:= 0 to n-1 do
    res[i] := round(fa[i].re + 0.5);
  carry:=0;
  for i:= 0 to n-1 do begin
    res[i]:=res[i]+carry;
    carry := round(res[i] / 10);
    res[i] :=res[i] div 10;
  end;   
  k:=n-1;
  while (res[k] = 0) do
    dec(k);
  setlength(res,k+1);  
end;