Системы линейных уравнений

Системы линейных уравнений Ax=b можно поделить на два класса. Первый класс - это системы с квадратной невырожденной матрицей A, имеющие единственное решение. Вторым классом являются системы с прямоугольной матрицей A произвольного размера, возможно вырожденной. Эти системы обычно решаются с минимизацией квадрата невязки. Пакет ALGLIB предлагает подпрограммы для поиска решений обеих проблем.

Содержание

    Невырожденные системы: матрица общего вида
    Невырожденные системы: симметричная/Эрмитова положительно определенная матрица
    Возможно вырожденные системы
    Быстродействие
    Итеративное улучшение
    Manual entries

Невырожденные системы: матрица общего вида

Для решения систем с квадратной невырожденной матрицей общего вида могут быть использованы следующие подпрограммы:

RMatrixSolve
RMatrixSolveM
RMatrixLUSolve
RMatrixLUSolveM
RMatrixMixedSolve
RMatrixMixedSolveM
CMatrixSolve
CMatrixSolveM
CMatrixLUSolve
CMatrixLUSolveM
CMatrixMixedSolve
CMatrixMixedSolveM

Подпрограммы могут быть поделены на следующие группы:

Решающие задачу с вещественной матрицей (имя начинается на R) или с комплексной матрицей (имя начинается на C).
Решающие задачу с одной правой частью или с несколькими правыми частями (имя заканчивается на M).
Принимающие только матрицу системы, только её LU разложение (имя содержит LU), либо и то, и то (имя содержит Mixed).

Последнее нуждается в уточнении:

Подпрограммы, принимающие матрицу системы, самостоятельно строят LU-разложение матрицы за время O(N³) и находят решение. К полученному решению может применяться итеративное улучшение (если параметр RFS установлен в True) с использованием арифметики повышенной точности (см. ниже).
Подпрограммы, принимающие LU-разложение матрицы системы, могут найти решение задачи за время O(N²), т.к. им не надо строить LU-разложение. Однако процедура итеративного улучшения в этом случае недоступна, т.к. итеративное улучшение требует точного знания оригинальной матрицы системы.
Если в подпрограмму передается и оригинальная матрица, и её LU-разложение, то мы можем сочетать преимущества обеих вариантов, приведенных выше: высокую точность первого варианта и быструю скорость второго варианта.

Невырожденные системы: симметричная/Эрмитова положительно определенная матрица

Для решения систем с симметричной/Эрмитовой положительно определенной матрицей общего вида могут быть использованы следующие подпрограммы:

SPDMatrixSolve
SPDMatrixSolveM
SPDMatrixCholeskySolve
SPDMatrixCholeskySolveM
HPDMatrixSolve
HPDMatrixSolveM
HPDMatrixCholeskySolve
HPDMatrixCholeskySolveM

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

Возможно вырожденные системы

Для решения систем с вещественной прямоугольной матрицей общего вида (возможно вырожденной) служит подпрограмма RMatrixSolveLS, основанная на сингулярном разложении матрицы системы. Если матрица система вырождена или близка к вырожденной, то минимизируется невязка. Подпрограмма поддерживает итеративное улучшение. Текущая версия ALGLIB не содержит версии этой подпрограммы для комплексных систем или для систем с множественными правыми частями.

RMatrixSolveLS solves systems with rectangular possibly degenerate matrix. Is system matrix is rank deficient, solution in a least squares sense is searched for. Subroutine supports iterative refinement. Current version of ALGLIB doesn't include complex version of this subroutine or version which accept multiple right-hand parts.

Быстродействие

На быстродействие линейного солвера влияет несколько факторов:

Быстродействие треугольной факторизации системы. Для вычисления треугольной факторизации используются оптимизированные алгоритмы, эффективно использующие все уровни иерархии кэша процессора. Вместе с тем, следует отметить, что эффективность алгоритмов сильно зависит от того, какую версию ALGLIB вы используете. Так, 100% .NET версия значительно уступает по своей эффективности версии на C++.
Количество правых частей. Солвер оптимизирован для систем с одной правой частью или же с малым количеством правых частей. Если же объем правой части сравним с объемом матрицы системы, и обе матрицы велики и не помещаются в кэш, то на быстродействие начинает влиять отсутствие блокирования для различных уровней кэша. В таких случаях вы можете самостоятельно осуществить факторизацию матрицы с использованием кэш-эффективного алгоритма из модуля trfac, после чего воспользоваться одной из TRSM-подпрограмм модуля ablas для кэш-эффективного решения полученной системы.
Итеративное улучшение. Для систем с большой матрицей и одной/несколькими правыми частями его стоимость незначительна в сравнении со стоимостью факторизации. Однако если размер матрицы мал, либо велико количество правых частей, то итеративное улучшение может составить значительную часть общего времени решения системы.
Копирование матрицы в ходе решения. Алгоритмы, принимающие LU-разложение матрицы, не требуют значительных количеств дополнительной памяти. Если же в алгоритм передается сама матрица системы, то требуется выделение дополнительной памяти для создания её копии. Для очень больших матриц, сравнимых по объему с доступной оперативной памятью, это может оказать существенное влияние на быстродействие.

Итеративное улучшение

Итеративное улучшение - это метод, позволяющий улучшить качество решения системы линейных уравнений, имея одновременно матрицу системы A и её треугольную факторизацию. После того, как система решена и получено решение x₀, мы вычисляем остаток r₀=b-Ax₀ и используем его для уточнения решения: Ad₀=r₀, x₁=x₀+d₀ Операция может повторяться несколько раз. В пакете ALGLIB используется двухпроходное улучшение.

Замечание #1
Итеративное улучшение позволяет уменьшить ошибки, возникающие из-за округления, и сделать решение настолько точным, насколько это позволяет обусловленность матрицы A. Однако качество полученного результата зависит от того, с какой точностью вычисляется b-Ax₀. Для того, чтобы итеративное улучшение оказалось стоящим усилий, требуется вычислять эту разность с существенно более высокой точностью, чем используемая при решении системы. В пакете ALGLIB для итеративного улучшения используется переносимый матричный мультипликатор, позволяющий вычислять b-Ax₀ с максимальной точностью, имеющей смысл при наличии в A и b ошибок величиной не более 1 ulp. "Переносимый" в данном контексте обозначает, что он не использует специфичных для той или иной архитектуры расширений (например, 80-битных расширений в архитектуре Intel) и не содержит низкоуровневых ассемблерных вставок.