LEC-28 (Материалы к лекциям), страница 2

Над каждым элементом а_ij выполняется самое большое 2n² арифметических операций, поэтому С=n² - самое большое значение, которое следует рассматривать. Вообще C=n достаточно надежный выбор.

Подведем итоги:

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

2. Не известны способы масштабирования матриц в общем случае, и существуют матрицы, которые не могут быть удовлетворительно масштабированы.

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

4. Для получения LU - разложения исключением Гаусса требуется порядка операций сложения и умножения.

Теперь обратимся к методу Холесского - симметричному варианту Гауссова исключения для симметричных положительно определенных матриц.

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

Ax=b, где А - NN симметричная положительно определенная матрица коэффициентов,

b - вектор длины N, называемый правой частью, а

х - вектор длины N, компоненты которого нужно вычислить.

Применения к А метода Холесского приводит к треугольному разложению

А=L L^T, где L - нижняя треугольная матрица с положительными диагональными элементами.

Верхний индекс Т указывает на операцию транспонирования.

Таким образом имеем: L L^T х=b

Замена y= L^Tх показывает, что х можно получить, решая треугольные системы:

Ly=b

L^T х=b

Поскольку L - треугольная матрица, то строка K может содержать ненулевые элементы только в первых K столбцах.

В частности, первая строка L содержит только один ненулевой элемент, так, что первая компонента y может быть легко вычислена. Вторая строка L содержит только два ненулевых элемента, следовательно, она связывает только две первые компоненты y; но так как первая компонента y известна, то теперь может быть вычислена и вторая. Таким образом, решаем систему уравнений относительно y.

Далее, поскольку матрица L^Т также треугольная, аналогичный процесс используется для отыскания х из второго уравнения, только теперь требуется обратный порядок операций, так как - L^Т верхняя, а не нижняя треугольная матрица.

Рассмотрим сначала случай 22 матрицы

А= L=

Потребуем, чтобы А=L L^T

L^T= : L L^T=

a₁₁=l₁₁²; a₁₂=l₁₁l₂₁; a₂₁=l₂₁l₁₁; a₂₂=l₂₁²+ l₂₂²

В силу положительной определенности матрицы А, элемент a₁₁>0 и извлечение квадратного корня l₁₁= возможно.

Второй элемент определяется из уравнения для a₂₁ (который равен a₁₂, т.к. А - симметричная матрица).

Наконец, l₂₂ определяется из равенства :

l₂₂ =

Алгоритм разложения Холесского для симметричной матрицы:

Цикл по столбцам А

For K=1 to n do

Вычисление очередного вектора-строки, кроме диагонального элемента.

For i=1 to K-1 do

a_Ki=( a_Ki - )/a_ii

Конец цикла по i.

Вычисление диагонального элемента

Конец цикла по К.

Теперь множитель L расположен в нижнем треугольнике матрицы А.

Этот алгоритм может быть сделан более устойчивым заменой оператора, содержащего квадратный корень, следующим оператором:

IF t >0 then a_KK= else a_KK =0

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

12.2. Проблема упорядочения

Вернёмся к примеру 2 решения системы уравнения методом Гаусса.

Обращает на себя внимание тот факт, что в процессе разложения матрица может претерпевать заполнение (на месте нулевых элементов матрицы появляются ненулевые элементы).

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

Множитель Холесского для матрицы коэффициентов системы выглядит так:

L =

Решая систему Ly=b, получаем: y=

Затем, решая систему L^TX=y, находим: X=

Этот пример иллюстрирует наиболее важный факт, относящийся к применению метода Холесского для разреженной матрицы А: матрица обычно претерпевает заполнение.

Это значит, что L имеет ненулевые элементы в позициях, где в нижней треугольной части А стояли нули.

Предположим теперь, что мы перенумеровали переменные в соответствии с правилом

Xi  , i =1,2,...,5

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

Должно быть ясно, что эта перенумерация переменных и переупорядочение уравнений равносильны симметричной перестановке строк и столбцов А, причем та же перестановка применяется к b.

Эту новую систему обозначим через

Применяя к ней , как и прежде, метод Холесского, мы разложим А в произведение , где (с точностью до 3-х значащих цифр)

Решая Ly=b и L^Tx=y, получим решение x , которое есть всего лишь переупорядоченная форма вектора х.

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

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

Существуют три серьезных причины, по которым заполнение следует считать нежелательным явлением:

1. Необходимо отводить память для хранения возникших
ненулевых элементов. При этом степень заполнения может ока­заться весьма значительной. Наиболее крайний пример, часто
цитируемый в литературе - это матрица с заполненными первой строкой,
первым столбцом и диагональю и остальными элементами, рав­ными нулю: если выполнить исключение в первом столбце, то
все остальные позиции матрицы подвергнутся заполнению. Если
матрица А симметрична и положительно определена, то размеры
заполнения и позиции новых ненулевых элементов могут быть
определены до того, как фактически начнется исключение. Таким
образом, становится возможным заранее выделить необходимый
дополнительный объем памяти. Однако если матрица А незнако-
определена, то это не удается сделать, и проблема распределения
памяти существенно усложняется.

2.Время, затрачиваемое при выполнении факторизации на
компьютере, быстро увеличивается с ростом заполнения, так как при­ходится выполнять гораздо большее количество арифметических
операций. Обычно объем памяти и число операций вместе определяют, поддается ли вообще данная задача решению прямыми
методами. Круг задач, к которым применимы прямые методы, тесно
связан, таким образом, с нашей способностью сводить заполнение
к минимуму.

3. Границы ошибок увеличиваются вместе с ростом заполне­ния. Многие исследователи отмечают ключевую роль того факта, что каждый элемент матрицы А лишь несколько раз претерпевает изменение в процессе исключения, и поэтому чрезмерного накопления вычислительных ошибок не происходит. С другой стороны, сдерживание роста ошибок должно рассматриваться как одна из целей, преследуемых при разработке алгоритмов, уменьшающих заполнение; этой точке зрения уделяется недостаточно внимания в литературе по разре­женным матрицам. Чрезмерное заполнение может повлечь за собой слишком сильное накопление вычислительных ошибок, что вынудит применить итерационный метод вместо прямого при решении системы линейных алгебраических уравнений большого размера.

К настоящему времени известно достаточно много об источ­никах возникновения заполнения и для его уменьшения разрабо­таны эффективные процедуры. Все они основаны на использовании свободы в выборе главных элементов, которые в случае положи­тельно определенной матрицы А могут выбираться на диагонали в любом порядке. Можно использовать другую интер­претацию этого подхода, согласно которой строки и столбцы матрицы А переупорядочиваются любым способом, сохраняющим симметрию, а затем производится исключение с выбором главных элементов на диагонали в естественном порядке. Хорошо известен тот факт, что размеры заполнения чрезвычайно существенно за­висят от того, какая именно выбрана перестановка.

Так, если рас­смотреть упомянутый выше «крайний» пример матрицы, все эле­менты которой равны нулю, кроме первой строки, первого столбца и диагонали, то при перестановке первого и последнего ее столб­цов и первой и последней строк ненулевыми в полученной матрице будут последняя строка, последний столбец и диагональ, так что при исключении вообще не возникает ни одного нового не­нулевого элемента.

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

Говорят, что упорядочение является оптимальным в смысле заполнения, если оно приводит к наименьшему возможному за­полнению. Упорядочение считается оптимальным в смысле затрат арифметических операций, если число операций, выполняемых при обработке переупорядоченной матрицы, минимально. Для матрицы А порядка п существует n! различных упорядочений, из которых одно или более оказывается оптимальным в смысле за­полнения, и несколько упорядочений (не менее одного) оказы­ваются оптимальными в смысле минимизации числа арифмети­ческих операций. Можно было бы пожелать найти такие опти­мальные упорядочения, однако, к сожалению, задача их отыска­ния представляется чрезвычайно трудной и о существовании эффек­тивных алгоритмов, пригодных для этой цели, не известно ни­чего. Существующие процедуры являются эвристическими и обычно основываются на попытках отыскания упорядочений, обеспечи­вающих понижение заполнения и числа арифметических операций, но не гарантирующих достижения точного минимума этих величин.