Решение систем линейных уравнений методом ортогонализации

Идея решения систем линейных уравнений методом ортогонализации состоит в следующим. Пусть даны матрица A (n*m) и вектор b (n), где n <= m. Надо найти такой x (m), чтобы Ax = b. Если n < m, то такая система называется недоопределённой ( к-во уравнений меньше, чем к-во переменных ). Перенесём вектор b в левую часть уравнения и будем рассматривать (m+1)-мерные строки. Тогда решением будет вектор x' у которого (m+1)-ая координата равна -1 и все скалярные произведения на строки полученной матрицы равны 0. Отсюда следует, что нам нужно найти в (m+1)-мерном пространстве вектор ортогональный каждому из n векторов. Такой вектор всегда существует. Если у него (m+1)-ая координата не равна 0, то умножая этот вектор на скаляр получим вектор x'.

Ортогональный вектор можно находить различными способами. Вот один из вариантов. Добавим к матрице строку у которой все координаты, кроме последней, равны 0 и применим к n+1 строкам процедуру ортогонализации Грама-Шмидта. В результате получим множество ортогональных векторов. Последний из них приводим к необходимому виду, если это возможно. Этот алгоритм находит всегда одно решение, даже в случае недоопределённой системы:

bool slu_ortho ( unsigned int n, unsigned int m, const double * const * a, const double * b, double * x );
Если же воспользоваться алгоритмом orthogonalizationH1 из раздела Ортогонализация векторов, то получим следующие две функции:
bool slu_orthoH1 ( unsigned int n, unsigned int m, const double * const * a, const double * b, double * x );

unsigned int slu_orthoH1 ( const Matrix & a, const double * b, Matrix & x );
Первая из них находит всегда одно решение, а вторая - в зависимости от ранга матрицы А и возвращает к-во найденных независимых решений.

Исходники находятся в файле mathem.cpp.

Наверх