В этом примере показано, как вычислить QR-разложение матриц с использованием аппаратного кода MATLAB ® в Simulink ®.
Чтобы решить систему уравнений или вычислить решение методом наименьших квадратов для матричного уравнения AX = B с использованием QR разложения, вычислите R и Q 'B, где QR = A, и RX = Q' B. R - верхняя треугольная матрица, а Q - ортогональная матрица. Если вы просто хотите Q и R, то установите B в качестве единичной матрицы.
В этом примере R вычисляется из матрицы A путем применения преобразований Гивенса с использованием алгоритма CORDIC. C = Q 'B вычисляется из матрицы B путем применения тех же преобразований Гивенса. Алгоритм использует только итеративные сдвиги и дополнения для выполнения этих вычислений.
Дополнительные сведения об алгоритме, используемом в этом примере, см. в разделе Выполнение факторизации QR с помощью CORDIC
Модель Simulink, используемая в этом примере:
fxpdemo_real_4x4_systolic_array_QR_model
Для ввода собственных входных матриц A и B откройте слева параметры соответствующих включенных блоков подсистемы. После моделирования модель возвращает вычисленные выходные матрицы C = Q 'B и R в рабочую область. Можно указать количество итераций CORDIC в параметрах блока подсистемы транспонирования (Q) * B, R 4x4 Real CORDIC систолического массива. Если входные данные являются фиксированными точками, то число итераций CORDIC должно быть меньше длины слова. Точность вычислений улучшает один бит для каждой итерации, вплоть до длины слова входных данных.
Эта модель будет работать с фиксированными, двойными и одиночными типами данных.
Чтобы увидеть, как алгоритм выполняет факторизацию, посмотрите под маской подсистемы транспонирования (Q) * B, R 4x4 Real CORDIC Sistolic Array. Аннотации указывают, какие строки матрицы обрабатываются для обнуления субдиагональных элементов. Систолический массив настраивается для матрицы A 4 на 4, но может быть расширен до любого размера, следуя тому же шаблону. Эта реализация работает только с реальными входными матрицами.
Чтобы увидеть код MATLAB в блоке MATLAB Function, который выполняет преобразования Givens с помощью CORDIC, продолжайте искать под масками блоков.
В этом примере число строк и столбцов A должно быть равно 4. Матрица B должна содержать 4 строки и любое количество столбцов.
Первый этап решения матричного уравнения AX = B состоит в вычислении RX = Q 'B, где R - верхнетреугольный, Q - ортогональный и Q * R = A.
Следующие входные данные являются типами с плавающей запятой двойной точности, поэтому установите число итераций равным 52, что является числом битов в мантиссе двойной.
format NumberOfCORDICIterations = 52; A = 2*rand(4,4)-1; B = 2*rand(4,4)-1;
Смоделировать модель для вычисления R и C = Q 'B.
sim fxpdemo_real_4x4_systolic_array_QR_model
R
C
R =
1.5149 -0.0519 1.7292 -0.3224
0 0.9593 -0.0259 -0.0879
0 0 0.2565 1.0888
0 0 0 -0.6429
C =
0.5942 -0.2382 0.0676 -0.9370
-0.8887 0.6146 -0.5758 0.3051
0.1725 0.7339 0.5409 0.5374
0.8540 1.1078 -0.2183 -0.5620
Убедитесь, что обратная замена на R и C = Q 'B дает те же результаты, что и MATLAB.
X = R\C X_should_be = A\B
X =
-7.1245 -12.1131 -0.6637 1.4236
-0.8779 0.7572 -0.5511 0.3545
6.3113 10.1759 0.6673 -1.6155
-1.3283 -1.7231 0.3396 0.8741
X_should_be =
-7.1245 -12.1131 -0.6637 1.4236
-0.8779 0.7572 -0.5511 0.3545
6.3113 10.1759 0.6673 -1.6155
-1.3283 -1.7231 0.3396 0.8741
Норма разницы между встроенным MATLAB и решением CORDIC QR должна быть небольшой.
norm(X - X_should_be)
ans = 3.6171e-14
Для вычисления Q и R установите B равным единичной матрице. Когда B равно единичной матрице, то Q = C '.
NumberOfCORDICIterations = 52; A = 2*rand(4,4)-1; B = eye(size(A,1),'like',A); sim fxpdemo_real_4x4_systolic_array_QR_model Q = C';
Теоретическое QR разложение QR = = A, поэтому разница между вычисленными QR и A должна быть небольшой.
norm(Q*R - A)
ans = 2.2861e-15
Декомпозиция QR уникальна только до знаков строк R и столбцов Q. Можно сделать уникальную декомпозицию QR, сделав все диагональные элементы R положительными.
D = diag(sign(diag(R))); Qunique = Q*D Runique = D*R
Qunique =
-0.3086 0.1224 -0.1033 -0.9376
-0.6277 -0.7636 -0.0952 0.1174
-0.5573 0.3930 0.7146 0.1559
0.4474 -0.4975 0.6852 -0.2877
Runique =
1.4459 -0.8090 0.1547 0.3977
0 1.1441 0.0809 -0.2494
0 0 0.8193 0.1894
0 0 0 0.4836
Затем можно сравнить вычисленный QR из модели с функцией builtin MATLAB qr.
[Q0,R0] = qr(A); D0 = diag(sign(diag(R0))); Q0 = Q0*D0 R0 = D0*R0
Q0 =
-0.3086 0.1224 -0.1033 -0.9376
-0.6277 -0.7636 -0.0952 0.1174
-0.5573 0.3930 0.7146 0.1559
0.4474 -0.4975 0.6852 -0.2877
R0 =
1.4459 -0.8090 0.1547 0.3977
0 1.1441 0.0809 -0.2494
0 0 0.8193 0.1894
0 0 0 0.4836
Используйте типы входных данных с фиксированной точкой для создания эффективного кода HDL для устройств ASIC и FPGA.
Дополнительные сведения о выборе типов данных с фиксированной точкой, которые не будут переполняться, см. в примере Выполнение факторизации QR с помощью CORDIC.
Можно выполнить множество тестовых входов через модель, создав 3-мерные массивы A и B.
n_test_inputs=100;
Следующий раздел определяет случайные входы для матриц A и B, которые масштабируются между -1 и 1, поэтому установите длину слова с фиксированной точкой в 18 бит и длину дроби в 14 бит, чтобы обеспечить рост факторизации QR и промежуточных вычислений в алгоритме CORDIC.
word_length = 18; fraction_length = 14;
Число итераций CORDIC с наилучшей точностью - это длина слова минус один. Если для числа итераций CORDIC установлено значение меньше word_length - 1, то задержка и такты следующего готового сигнала будут короче, но будут менее точными. Количество итераций CORDIC не должно быть больше, поскольку созданный код не поддерживает сдвиги, превышающие длину слова типа с фиксированной точкой.
NumberOfCORDICIterations = word_length - 1
NumberOfCORDICIterations =
17
Случайные тестовые входы объединяются так, что в момент времени k входами являются A (:,:, k) и B (:,:, k). Каждый элемент A и B является однородной случайной величиной между -1 и + 1.
A = 2*rand(4,4,n_test_inputs)-1;
Выберите B в качестве единичной матрицы, так что Q = C '.
B = eye(4); B = repmat(B,1,1,n_test_inputs);
Отливка А в фиксированной точке и отливка В, как А.
A = fi(A,1,word_length,fraction_length);
B = cast(B,'like',A);
sim fxpdemo_real_4x4_systolic_array_QR_model
norm_error = zeros(1,size(R,3)); for k = 1:size(R,3) Q_times_R_minus_A = double(C(:,:,k))'*double(R(:,:,k)) - double(A(:,:,k)); norm_error(k) = norm(Q_times_R_minus_A); end
Ошибки должны быть порядка 10 ^ -3.
clf plot(norm_error,'o-') grid on title('norm(QR - A)')
%#ok<*NASGU,*NOPTS>