Реализуйте эффективный оборудованием комплексный частично-систолический Q-less QR с упущением фактора

Открыть скрипт

В этом примере показано, как использовать эффективное оборудованием Комплексное Частично-систолическое Разложение Q-less QR с блоком Forgetting Factor.

Разложение Q-less QR с упущением фактора

Комплексное Частично-систолическое Разложение Q-less QR с блоком Forgetting Factor реализует следующую рекурсию, чтобы вычислить верхний треугольный множитель R непрерывной потоковой передачи n-1 векторы-строки (k, :) с помощью упущения фактора $\alpha$ . Это - как будто матрица А бесконечно высока. Упущение включает область значений $0<\alpha<1$ , мешает ему объединяться без связанного.

$
\begin{array}{ccl}
R_0 & =& \mbox{zeros}(n,n)\\[2ex]
[\sim\;,\;R_1] &=& \mbox{qr}\left(\alpha\left[\begin{array}{c}R_0\\
A(1,:)\end{array}\right],\; 0\right)\\[4ex]
[\sim\;,\;R_2] &=& \mbox{qr}\left(\alpha\left[\begin{array}{c}R_1\\
A(2,:)\end{array}\right],\; 0\right)\\[4ex]
\vdots\\[4ex]
[\sim\;,\;R_k] &=& \mbox{qr}\left(\alpha\left[\begin{array}{c}R_{k-1}\\
A(k,:)\end{array}\right],\; 0\right)\\[4ex]
\vdots\\[4ex]
\end{array}
$

Задайте системные параметры

n длина векторов-строк (k, :) и количество строк и столбцов в R.

n = 5;

m эффективные количества строк, чтобы объединяться.

m = 100;

Используйте fixed.forgettingFactor функция, чтобы вычислить фактор упущения в зависимости от количества строк, которые вы интегрируете.

forgettingFactor = fixed.forgettingFactor(m)

forgettingFactor =

    0.9950

precisionBits задает количество битов точности, требуемой для Разложения QR. Установите это значение согласно системным требованиям.

precisionBits = 24;

В этом примере матрица А с комплексным знаком создается таким образом, что величина действительных и мнимых частей ее элементов меньше чем или равна одному, таким образом, максимальное возможное абсолютное значение любого элемента $|1+1i|=\sqrt{2}$ . Ваши собственные системные требования зададут, каковы те значения. Если вы не знаете то, что они, и A является входом фиксированной точки к системе, то можно использовать upperbound функция, чтобы определить верхние границы фиксированных точек A.

max_abs_A верхняя граница на максимальном элементе массива величины.

max_abs_A = sqrt(2);

Выберите Fixed-Point Types

Используйте fixed.qlessqrFixedpointTypes функция, чтобы вычислить фиксированные точки.

T = fixed.qlessqrFixedpointTypes(m,max_abs_A,precisionBits)

T = 

  struct with fields:

    A: [0x0 embedded.fi]

T.A фиксированная точка, вычисленная для преобразования к R, оперативному так, чтобы это не переполнялось.

T.A

ans = 

[]

          DataTypeMode: Fixed-point: binary point scaling
            Signedness: Signed
            WordLength: 30
        FractionLength: 24

Задайте параметры симуляции

Создайте случайную матрицу А, чтобы содержать конкретное количество входных параметров.

numInputs количество входных строк A (k, :) для этого примера.

numInputs = 500;
rng('default')
A = fixed.example.complexUniformRandomArray(-1,1,numInputs,n);

Бросьте входные параметры к типам, определенным fixed.qlessqrFixedpointTypes.

A = cast(A,'like',T.A);

Бросьте фактор упущения к фиксированной точке с тем же размером слова как A и масштабирование лучшей точности.

forgettingFactor = fi(forgettingFactor,1,T.A.WordLength);

Выберите время остановки для симуляции, которая достаточно долга, чтобы обработать все входные параметры от A.

stopTime = 4*numInputs*T.A.WordLength;

Откройте модель

model = 'ComplexPartialSystolicQlessQRForgettingFactorModel';
open_system(model);

Установите переменные в рабочем пространстве модели

Используйте функцию помощника setModelWorkspace добавить переменные, заданные выше к рабочему пространству модели.

fixed.example.setModelWorkspace(model,'A',A,'n',n,...
    'forgettingFactor',forgettingFactor,...
    'stopTime',stopTime);

Симулируйте модель

out = sim(model);

Проверьте точность Выхода

Задайте матрицу $A_k$ можно следующим образом

$A_k = \left[\begin{array}{cccc}\alpha^k \\& \alpha^{k-1} \\
& & \ddots\\&
& & \alpha\end{array}\right]A(1:k,\; :).$

Затем с помощью формулы для расчета $k$ th выход $R_k$ и то, что $[Q,R]=\mbox{qr}(A,0) \Rightarrow A'A = R'Q'QR = R'R$ , можно показать это

$$A_k'A_k = R_k'R_k.$$

Таким образом проверить выход, различие между $A_k'A_k$ и $R_k'R_k$ должно быть малым.

Выберите последний выход симуляции.

R = double(out.R(:,:,end))

R =

  Columns 1 through 4

   7.8315 + 0.0000i   0.2537 + 0.1037i  -0.1812 + 0.1835i   0.7054 + 0.2411i
   0.0000 + 0.0000i   7.7917 + 0.0000i  -0.2939 - 0.6734i  -0.2813 - 0.2824i
   0.0000 + 0.0000i   0.0000 + 0.0000i   8.2110 + 0.0000i  -0.2445 - 0.0531i
   0.0000 + 0.0000i   0.0000 + 0.0000i   0.0000 + 0.0000i   8.0730 + 0.0000i
   0.0000 + 0.0000i   0.0000 + 0.0000i   0.0000 + 0.0000i   0.0000 + 0.0000i

  Column 5

  -0.2351 - 0.1305i
  -0.4660 - 0.6805i
  -0.2956 + 0.5699i
   0.2053 - 0.0133i
   8.0363 + 0.0000i

Проверьте, что R верхний треугольный.

isequal(R,triu(R))

ans =

  logical

   1

Проверьте, что диагональ больше или равна нулю.

diag(R)

Синхронизируйте последний выход R с входом путем нахождения количества входных параметров, которые произвели его.

A = double(A);
alpha = double(forgettingFactor);
relative_errors = nan(1,n);
for k = 1:numInputs
    A_k = alpha.^(k:-1:1)' .* A(1:k,:);
    relative_errors(k) = norm(A_k'*A_k - R'*R)/norm(A_k'*A_k);
end

k количество входных параметров (k, :) который произвел последний R.

k = find(relative_errors==min(relative_errors),1,'last')

k =

   497

VerifyThat

$$A_k'A_k = R_k'R_k$$

с небольшой относительной погрешностью.

A_k =  alpha.^(k:-1:1)' .* A(1:k,:);
relative_error = norm(A_k'*A_k - R'*R)/norm(A_k'*A_k)

relative_error =

   6.3862e-06

Подавите mlint предупреждения в этом файле.

%#ok<*NOPTS>

Документация