Когда ортогональные преобразования превосходят исключение Гаусса?

22

Как мы знаем, методы ортогональных преобразований (повороты Гивенса и отражения Хаусхолдера) для систем линейных уравнений более дороги, чем устранение по Гауссу, но теоретически обладают более хорошими свойствами устойчивости в том смысле, что они не изменяют число условий системы. Хотя я знаю только один академический пример матрицы, которая испорчена исключением Гаусса с частичным поворотом. И есть общее мнение, что на практике такого рода поведение вряд ли встретится (см. Примечания к этой лекции [pdf] ).

Итак, где мы будем искать ответ по теме? Параллельные реализации? Обновление? ..

faleichik
источник

Ответы:

24

точность

Трефетен и Шрайбер написали отличную статью « Стабильность по Гауссу в среднем случае» , в которой обсуждается вопрос о точности вашего вопроса. Вот несколько его выводов:

  1. «Для QR разложения с или без колонки поворота, средний максимальный элемент остаточной матрицы , тогда как для исключения Гаусса это O ( п ) . Это сравнение показывает , что исключение Гаусса умеренно неустойчивым, но нестабильность может быть обнаружена только для очень больших матричных задач, решаемых с низкой точностью. Для большинства практических задач элиминация по Гауссу в среднем очень стабильна. "(Выделение мое)O(n1/2)O(n)

  2. «После первых нескольких этапов исключения из Гаусса остальные матричные элементы примерно нормально распределены, независимо от того, были ли они начаты таким образом».

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

Спектакль

Для квадратных вещественных матриц, LU с частичной поворота требуется примерно провалов, тогда как Хаусхолдера на основе QR - требуется примерно 4 / 3 л 32/3n34/3n3 провалов. Таким образом, для достаточно больших квадратных матриц QR-факторизация будет в два раза дороже, чем LU-факторизация.

Для матриц, где т п , LU с частичным поворотом требует м п 2 - п 3 / 3 флопа, по сравнению с QR -- х 2 м н 2 - 2 л 3 / 3 (который по - прежнему вдвое больше , чем LU факторизации). Тем не менее , неожиданно часто приложения создают очень высокие узкие матрицы ( m nm×nmnmn2n3/32mn22n3/3mn ), и Demmel et al. приятной работы, избегая общения параллельной и последовательной QR-факторизации, который (в разделе 4) обсуждает умный алгоритм, который требует, чтобы сообщения отправлялись только при использовании p процессоров, в отличие от n log p сообщений традиционных подходов. Суть в том, что O ( n 3 log p ) выполняются дополнительные провалы, но для очень малых n это часто предпочтительнее, чем стоимость задержки отправки большего количества сообщений (по крайней мере, когда необходимо выполнить только одну факторизацию QR).logppnlogpO(n3logp)n

Джек Полсон
источник
10

Я удивлен, что никто не упомянул линейные задачи наименьших квадратов , которые часто встречаются в научных вычислениях. Если вы хотите использовать метод исключения Гаусса, вы должны сформировать и решить нормальные уравнения, которые выглядят так:

ATAx=ATb,

где - матрица точек данных, соответствующих наблюдениям независимых переменных, x - вектор параметров, которые нужно найти, и bAxb - вектор точек данных, соответствующих наблюдениям зависимой переменной.

Как часто отмечает Джек Полсон, число условий является квадратом числа условий A , поэтому нормальные уравнения могут быть катастрофически плохо обусловлены. В таких случаях, хотя QR-и SVD-подходы медленнее, они дают гораздо более точные результаты.ATAA

Джефф Оксберри
источник
2
n3AHA2/3n36n3
1
В дополнение к стабильности, гарантируемой использованием ортогональных преобразований, большое преимущество SVD состоит в том, что декомпозиция обеспечивает свою собственную проверку условий, поскольку отношение наибольшего к наименьшему единственному значению является точно (2-нормальным) числом условия. Для других разложений использование оценщика условий (например, Хагера-Хайама), хотя и не так дорого, как собственно разложение, несколько «привязано».
JM
1
@JackPoulson Просто из любопытства, у вас есть ссылка на счет флопа для SVD? Из того, что я могу сказать из краткого обзора в Golub & Van Loan (стр. 254, 3-е издание), константа может показаться выше для использования SVD при решении задач наименьших квадратов, но я могу ошибаться. Заранее спасибо.
OscarB
1
8/3n3A=FBGHCB=UΣVHx:=(G(V(inv(Σ)(UH(FHb)))))O(n2)CO(n2)
1
σ1σn
3

Как вы измеряете производительность? Скорость? Точность? Стабильность? Быстрый тест в Matlab дает следующее:

>> N = 100;
>> A = randn(N); b = randn(N,1);
>> tic, for k=1:10000, [L,U,p] = lu(A,'vector'); x = U\(L\b(p)); end; norm(A*x-b), toc
ans =
   1.4303e-13
Elapsed time is 2.232487 seconds.
>> tic, for k=1:10000, [Q,R] = qr(A); x = R\(Q'*b); end; norm(A*x-b), toc             
ans =
   5.0311e-14
Elapsed time is 7.563242 seconds.

Таким образом, решение одной системы с LU-декомпозицией примерно в три раза быстрее, чем ее решение с помощью QR-декомпозиции, за счет половины точности до десятичного знака (этот пример!).

Pedro
источник
Любые из предложенных вами достоинств приветствуются.
фалейчик
3

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

Рассмотрим пространство всех матриц. Эти методы прекрасно работают практически везде. То есть 99,999 ...% всех матриц, которые можно создать, не будут иметь проблем с нестабильными методами. Существует только очень небольшая часть матриц, для которых GE и другие будут испытывать трудности.

Проблемы, которые волнуют исследователей, как правило, в этой небольшой части.

Мы не строим матрицы случайно. Мы строим матрицы с очень особыми свойствами, которые соответствуют совершенно особым неслучайным системам. Эти матрицы часто плохо обусловлены.

Геометрически вы можете рассмотреть линейное пространство всех матриц. Существует подпространство нулевого объема / меры сингулярных матриц, прорезающих это пространство. Многие проблемы, которые мы строим, сгруппированы вокруг этого подпространства. Они не распределены случайно.

В качестве примера рассмотрим уравнение теплопроводности или дисперсию. Эти системы имеют тенденцию удалять информацию из системы (все начальные состояния тяготеют к одному конечному состоянию), и в результате матрицы, которые описывают эти уравнения, чрезвычайно необычны. Этот процесс очень маловероятен в случайной ситуации, но распространен в физических системах.

MRocklin
источник
2
Если линейная система изначально плохо обусловлена, то независимо от того, какой метод вы используете: декомпозиция LU и QR даст неверные результаты. QR может победить только в тех случаях, когда процесс исключения из Гаусса «портит» хорошую матрицу. Основная проблема заключается в том, что практические случаи такого поведения не известны.
фалейчик
Для большинства научных приложений мы обычно получаем матрицы, которые являются разреженными, симметричными, положительно определенными и / или диагонально доминирующими. За очень немногими исключениями, в матрице есть структура, которая позволяет нам использовать определенные методы по сравнению с традиционным устранением по Гауссу.
Павел
@Paul: С другой стороны, плотное исключение Гаусса - то, где большая часть времени проводится в мультифронтальном методе для разреженных несимметричных матриц.
Джек Полсон,
6
@Paul Это просто неправда, что «большинство приложений производят SPD / диагонально доминирующие матрицы». Да, обычно есть какая-то эксплуатируемая структура, но несимметричные и неопределенные проблемы встречаются крайне часто.
Джед Браун
4
«За пятьдесят лет вычислений не было обнаружено никаких проблем с матрицами, которые возбуждают взрывную нестабильность, в естественных условиях». - Л. Н. Трефетен и Д. Бау. В своей книге они дают интересный вероятностный анализ.
JM