Обращение матриц на ПЛИС Опции

[jojo]

Обращение матриц на ПЛИС

Тут надо комплексную матрицу 16x16 обратить (эрмитову). Формат исходных данных int32, результат - не хуже float, лучше - double. Есть ли готовые мегафункции, не обязательно нелицензионные?

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

В нынешней реализации обращения на DSP используются квадратный корень и деление (вещественные), сложение, вычитание и умножение (комплексные).

Хочу уменьшить нагрузку на DSP и разместить обращение матрицы в относитильно защищенной от взлома ПЛИС.

[jojo]

#include <math.h>

Обращение выполняется за 4 этапа. С двумя последними мне все понятно, там "* + -" с плавающей точкой. И с регуляризацией понятно. Но во что обойдется реализация этапа 2 в ПЛИС. Как я понимаю, придется сделать значительный кусок АЛУ DSP, который делает sqrt и 1/x во флоате.


Вообще сколько стоит денег в Москве такое реализовать? Может, руководство согласится и вместо TS203 мы поставим ПЛИС и простенький DSP..

Но не люблю я эти Плисы. Тут есть много проектов, которые накрылись из-за того, что аффторы не справились из-за Плис.

A - массив 16x16 структур {float Re, float Im}.
d - число для регуляризации матрицы
N - 16, N2 - 32 (так автор алгоритма захотел)
S и K - временные массивы 16x16 структур {float Re, float Im}.

void inv(float *A, int N, int N2, float *S, float *K, float d)
{
int p,q,k,pr,pi,qr,qi,kr,ki;
float a,ar,ai;

//(Этап 0 приведен полностью) регуляризация
q=0;
for (p=0;p<N;p++)
{
A[p*N2+q]+=d;
q+=2;
}

//(Этап 1 приведен полностью) вычисление S
//p=1
a=A[0];
S[0]=sqrt(a);
S[1]=0.;
K[0]=1./S[0];
K[1]=0.;
//re и im
for (p=0;p<N2;p++)
{
S[p]=A[p]*K[0];
}

//p>1
pr=0;
pi=1;
for (p=1;p<N;p++)
{
pr+=2;
pi+=2;
a=A[p*N2+pr];
for (q=0;q<=p-1;q++)
{
a-=(S[q*N2+pr]*S[q*N2+pr]+S[q*N2+pi]*S[q*N2+pi]);
}


S[p*N2+pr]=(float)sqrt(a);
S[p*N2+pi]=0.;
K[p*N2+pr]=1./S[p*N2+pr];
K[p*N2+pi]=0.;


qr=pr;
qi=pi;
for (q=p+1;q<N;q++)
{
qr+=2;
qi+=2;
ar=A[p*N2+qr];
ai=A[p*N2+qi];
for (k=0;k<=p-1;k++)
{
ar-=(S[k*N2+pr]*S[k*N2+qr]+S[k*N2+pi]*S[k*N2+qi]);
ai-=(S[k*N2+pr]*S[k*N2+qi]-S[k*N2+pi]*S[k*N2+qr]);
}
S[p*N2+qr]=ar*K[p*N2+pr];
S[p*N2+qi]=ai*K[p*N2+pr];
}
} //(Этап 1 завершен)


//Этап 2 пропущен. Это перемножение матриц S и K
ar-=(S[p*N2+kr]*K[k*N2+qr]-S[p*N2+ki]*K[k*N2+qi]);
ai-=(S[p*N2+kr]*K[k*N2+qi]+S[p*N2+ki]*K[k*N2+qr]);
K[p*N2+qr]=ar*K[p*N2+pr];
K[p*N2+qi]=ai*K[p*N2+pr];

//Этап 3 пропущен. Это перемножение матриц K и транспонированной K
ar+=(K[p*N2+kr]*K[q*N2+kr]+K[p*N2+ki]*K[q*N2+ki]);
ai+=(K[p*N2+ki]*K[q*N2+kr]-K[p*N2+kr]*K[q*N2+ki]);
A[]=ar;
A[]=ai;
A[]=ar;
A[]=-ai;

return;
};