Помощь - Поиск - Пользователи - Календарь
Полная версия этой страницы: Быстрые алгоритмы нахождения обратной матрицы
Форум разработчиков электроники ELECTRONIX.ru > Cистемный уровень проектирования > Математика и Физика
Grizzzly
Подскажите, пожалуйста, самый быстрый алгоритм для нахождения обратной матрицы небольшого размера (6x6), разреженностью или симметрией не обладает.

Спасибо.

P.S. Или способ решения переопределенных линейных систем уравнений (8 уравнений, 6 неизвестных), обладающий наименьшей вычислительной сложностью.
thermit
Цитата(Grizzzly @ May 16 2016, 23:59) *
Подскажите, пожалуйста, самый быстрый алгоритм для нахождения обратной матрицы небольшого размера (6x6), разреженностью или симметрией не обладает.

Спасибо.

P.S. Или способ решения переопределенных линейных систем уравнений (8 уравнений, 6 неизвестных), обладающий наименьшей вычислительной сложностью.


Гаусс с выбором главного элемента. Потом идет хаусхолдер. Потом гивенс. Двое последних - численно устойчивые.
Grizzzly
Цитата(thermit @ May 17 2016, 01:11) *
Гаусс с выбором главного элемента. Потом идет хаусхолдер. Потом гивенс. Двое последних - численно устойчивые.

Спасибо.
Maverick
Цитата(Grizzzly @ May 17 2016, 08:23) *

посмотрите это
+
это
Grizzzly
Цитата(Maverick @ May 17 2016, 12:46) *

Супер!
Maverick
Цитата(Grizzzly @ May 17 2016, 12:48) *

для маленьких матриц и это
может пригодиться (по первой ссылке для 3х3), думаю можно увеличить до матрицы 6х6
Grizzzly
Цитата(Maverick @ May 17 2016, 13:10) *
думаю можно увеличить до матрицы 6х6

Спасибо! Попробую.
iiv
Цитата(Grizzzly @ May 17 2016, 02:59) *
P.S. Или способ решения переопределенных линейных систем уравнений (8 уравнений, 6 неизвестных), обладающий наименьшей вычислительной сложностью.

как я понимаю, исходная задача именно эта, а 6х6 получилась решением наименьших квадратов.

Если таки 8х6 решать, то надо сделать для этой матрицы QR разложение с выбором ведущего столбца, можно сделать даже Грамм-Шмидтом, тогда алгоритм будет помещаться в 5-6 строк. Если на в R на каком-то шаге на диагонали получаемое число будет меньше чем первое в некоторое эпсилон, которое у Вас больше ошибки входных данных, то Вам надобно регуляризовать, так как иначе получится на выходе лажа.

Самым кондовым методом тогда будет оборвать QR как только на диагонали будет такое, и получить Q R_r, где Q матрица 8хr, а R_r матрица rх6 размеров. Для транспонированной R_r повторить разложение, но до конца, получив Q P_2 R Q_r P_1, где P - матрицы перестановки. Если есть возможность не писать, а попользовать лапак, то просто вызвать DGESVD из оного и занулить не нужные сингулярные значения и решить именно так, но, лапак не на всех платформах компилится, а то, что я предложил, при необходимости в атмегу даже засунуть можно, а вот DGESVD уже точно не засунешь sm.gif

Если Вы до этого никогда не писали софт по линейной алгебре, то самопально написанный аогоритм Грамма-Шмидта или Хаусхолдера с перестановками и всем тем добром, что я написал, точно проиграет по скорости лапаку, хотя арифметическая сложность сингулярного разложения существенно выше Грамма-Шмидта или Хаусхолдера, это утверждение много раз мной было проверенно на куче как физтеховских, так и немецких студентах sm.gif
thermit
Цитата
iiv:
Если Вы до этого никогда не писали софт по линейной алгебре, то самопально написанный аогоритм Грамма-Шмидта или Хаусхолдера с перестановками и всем тем добром, что я написал, точно проиграет по скорости лапаку


Если только в лапаке магия какая реализована. Кроме того, кому нужно сингулярное разложение? Достаточно привести исходную систему к треугольной гауссом, хаусхолдером или гивенсом. После обратной подстановки получим решение. Или вычислить матрицу грама и решить полученную систему тем же методом. Или разложением холецкого, которое подкупает своей простотой, но содержит несколько подводных граблей.
iiv
Цитата(thermit @ May 18 2016, 03:48) *
Если только в лапаке магия какая реализована. Кроме того, кому нужно сингулярное разложение?.



ТС об обусловленности по входным данным ни слова не сказал, то есть скорей всего или не знает какая она, или вообще не знает что это. Поэтому и Гаусс, и Хаусхолдер с Гивенсом получат ему решение, но оно, скорей всего, будет далеким от того, что он ожидает и будет он плясать с бубном, пока годовой курс нумерики не пройдет. С SVD же он получит результат за один вызов.

PS: магия лапака в том, что у него все конвейерно и кеш-оптимизировано, а человек с этим до этого не сталкивавшийся быстро начать в этих терминах программировать не сможет, что я много раз воочию наблюдал
thermit
Цитата(iiv @ May 18 2016, 12:25) *
ТС об обусловленности по входным данным ни слова не сказал, то есть скорей всего или не знает какая она, или вообще не знает что это. Поэтому и Гаусс, и Хаусхолдер с Гивенсом получат ему решение, но оно, скорей всего, будет далеким от того, что он ожидает и будет он плясать с бубном, пока годовой курс нумерики не пройдет. С SVD же он получит результат за один вызов.

PS: магия лапака в том, что у него все конвейерной и кеш-оптимизировано, а человек с этим до этого не сталкивавшийся быстро начать в этих терминах программировать не сможет, что я много раз воочию наблюдал



Да нормальное решение выдадут хаусхолдер с гивенсом. Даже в случае совсем плохо обусловленной задачи. В крайнем случае в недрах алгоритмов получится деление на 0.
И пространства для плясок с бубном честно говоря здесь не вижу. Кроме того, кто сказал что это должно работать на пэцэ? Может человек все это в с66xx запихнуть хочет. А лапак, он такой лапак...
Swup
Вот интересная статья. Используется некоторая оптимизация поворотов Гивенса с реализацией на систолическом массиве.
https://www.researchgate.net/publication/22..._MIMO_Receivers

Если использовать ПЛИС, ну или что-то с большой параллельностью, то это очень быстрая реализация. Причем как раз эффективно для маленьких матриц.
Grizzzly
Цитата(Swup @ May 20 2016, 11:06) *
Вот интересная статья.

Огромнейшее спасибо! То что нужно. Еще и с таблицей количества операций. На работе тоже говорили про повороты Гивенса, нашел несколько русскоязычных статей, но эта в разы лучше.
serjj
Что-то не совсем понял исходную задачу. Ваша система имеет вид Ax = y, где размерность x - [6x1], y - [8x1], A - [8x6], задача оценить x, имея A и y?
Линейное решение получается через псевдоинверсию: x=(ATA)-1ATy, если A действительная и x=(AHA)-1AHy - если комплексная. В первом случае ATA симметричная матрица, во втором - AHA эрмитова. Эти матрицы обладают известной симметрией, поэтому на них разумно натравить алгоритм Холецкого, по скорости он выигрывает у методов, которые не используют специальные свойства матриц.
Grizzzly
Цитата(serjj @ Jun 18 2016, 00:19) *
Что-то не совсем понял исходную задачу. Ваша система имеет вид Ax = y, где размерность x - [6x1], y - [8x1], A - [8x6], задача оценить x, имея A и y?

Да, именно так. И решал через псевдообратную матрицу точно так же. С самим решением понятно. А вот с выбором алгоритма, отвечающего минимальной вычислительной сложности, были вопросы. Понял. Спасибо.

P.S. Задача свелась к чисто вещественной.
serjj
Цитата
P.S. Задача свелась к чисто вещественной.


Любопытно. Вы перешли от размерности N (6) к 2N (12) с последующим решением в поле вещественных чисел? Или у вас исходно только вещественный сигнал / матрица?
Grizzzly
Цитата(serjj @ Jun 19 2016, 13:30) *
Любопытно. Вы перешли от размерности N (6) к 2N (12) с последующим решением в поле вещественных чисел? Или у вас исходно только вещественный сигнал / матрица?

Изначально было 6 веещственных неизвестных и размерность матриц 4 (комплексные величины), перешел к 8 вещественным, получив переопределенную систему.
Для просмотра полной версии этой страницы, пожалуйста, пройдите по ссылке.
Invision Power Board © 2001-2024 Invision Power Services, Inc.