Численное решение задач метода наименьших квадратов - Лоуcон Ч.
Скачать (прямая ссылка):
4?)-0-«"'tT):
г
г.
41
ГЛАВА 10
ВЫЧИСЛЕНИЯ, ИСПОЛЬЗУЮЩИЕ ЭЛЕМЕНТАРНЫЕ ОРТОГОНАЛЬНЫЕ ПРЕОБРАЗОВАНИЯ
Мы переходим теперь к описанию численных алгоритмов для ортогонального разложения т X п-матрицы.4 (указанного в теореме 3.19) и решения задачи НК.
Так как преобразования Хаусхолдера и Гивенса имеют очень много приложений, мы сочли, что будет удобно и полезно с точки зрения унификации рассмотреть каждое из этих преобразований как самостоятельный объект. Оба преобразования будут сейчас разобраны в деталях. Эти модули в дальнейшем используются в описании других, более сложных вычислительных процедур.
Вычисление преобразования Хаусхолдера можно разбить на два этапа: 1) построение преобразования; 2) его применение к другим векторам.
Ортогональное т X т -преобразование Хаусхолдера можно представить в виде
Q = 1т + Ь-1ииТ, (10.1)
где и - m-вектор такой, что ||и|| Ф 0, и b - — \\и\\212.
В гл. 3 было показано, что для данного вектора и можно определить вектор и таким образом, чтобы выполнялось (3.2). На практике встречаются ситуации, когда Q нужно определить из условия
Чр-1
Qv =
/ m \ш o^v?+ Ъ vfj
Vp+l Vl-l
0
= У-
Построить такую матрицу Q можно было бы как произведение перестановки строк, обычного преобразования Хаусхолдера, описанного в лемме 3.1, и обратной перестановки строк. Мы находим, однако, более удобным рассматривать равенство (10.2) как определяющее цель нашего основного вычислительного модуля. Действие матрицы Q при преобразовании v в у можно описать с помощью трех неотрицательных целых параметров р, I, m следующим образом:
1) если р > 1, то компоненты с номерами 1..... р - 1 не должны
меняться;
42
2) компонента с номером р может измениться. Она называется главным элементом;
3) если р < I - 1, то компоненты с номерами р-н 1,/ - 1 не должны меняться;
4) если / < тп, то компоненты с номерами /, тп должны быть аннулированы.
Заметим, что должны выполняться соотношения
1 < р < тп, (10.3)
р < I. . (Ю.4)
Шаги численного процесса, который приводит к ортогОВяльной m X тп -матрице Q, отвечающей целым параметрам р, /, тп, МГОЮО свести в следующую схему:
m
s = -o(u2 + ? v2)112,
(10.5)
t = j + 1, vp > 0,
o = 1 ^ ¦ « -1, vp < 0,
m =o, / = 1.....p-i, (Ю.6)
up = vp-s, , . ' (10.7)
щ =0, J=p + 1,1, (10.8)
и, =v,. / = /.....m, (10.9)
ft =SUp, (10.10)
j-'um7", ft*0, (10.11)
= j/m + ft-'uH',
I An »
ft =0.
Тот факт, что матрица Q, определенная в (10.11), имеет нужные свойства, устанавливается следующими тремя леммами.
Лемма 10.12. Определяемые формулами (10.5) - (10.10) тп-век-тор и и скаляр b удовлетворяют равенству
ft=-||«||I/2. (10.13)
Доказательство. Имеем
Ни II2 = (vp - О2 + ? о2 =
/ = / '
= v2 - 2vps + s2 + ? и? •-i*(wp - j)--2uL --2ft.
j=i
Лемма доказана.
Лемма 10.14. Матрица Q, определяемом формулой (10.11),оргого-нальная.
43
Д о к а з а т е л ь с т в о. Проверка усяомю Q?Q я/1т проводится "в лоб" с использованием (10.11) и (10.13). /
Лемма 10.15. Пусть у- Qv. Тогда I
yi = и,, i = l, 1, f (10.16)
УР = s, (10.17)
У1 = «I. /=р + 1,...,/- 1, 1!" (Ю.18)
.у, - 0, i = l,...,m. (10.19)
Доказательство. Если v = 0, то лемма очевидным образом верна. При v Ф 0 легко проверяется с использованием (10.6) — (10.9), что uTv = —b. Но тогда вектор
у = Qv = v - и
удовлетворяет равенствам (10.16) - (10.19).
В реальных вычислениях ненулевые компоненты векторов и и у обычно размещают в массиве памяти, который прежде занимал вектор и. Исключением является р-я компонента массива; приходится выбирать, поместить ли в нее ир или ур. Мы остановимся на втором варианте, а для хранения ир заведем дополнительную ячейку. Величина b может быть вычислена в соответствии с (10.10) всякий раз. как она нужна.
После того как преобразование построено, обычно требуется применить его к некоторому набору m-векторов С\,cv, т.е. нужно вычислить векторы Cj = Qcj, j = 1,v. Используя определение Q в (10.11), это вычисление можно выполнять так:
t, = Ь-1(итс,), (10.20)
cj^ci+t/u, /=1,...,у. (10.21)
Все вышеизложенное будет теперь переформулировано в алгоритмической форме, пригодной для реализации машинной программой. Мы опишем алгоритм HI (р, /, пг, v, h. С, v) для построения и (по желанию пользователя) применения преобразования Хаусхолдера и алгоритм Н2(р,/,т, о, h,C,v) для применения (по желанию пользователя) построенного ранее преобразования Хаусхолдера.
Входной информацией алгоритма HI являются целые числа р,1,.тп и v, m-вектор v и при v > 0 массив С, содержащий m-векторы с/, j = l,...,i>.
Массив С может иметь либо размеры m X v, и в этом случае векторы с,-являются его столбцами, либо размеры v X m, и тогда Cf будут его строками. В описании алгоритмов HI и Н2 мы не будем различать зги два возможных способа хранения. Однако при ссылках на эти алгоритмы в последующем тексте книги хранение по столбцам будет считаться стандартным методом хранения. В тех же случаях, когда векторы с,-, к которым применяется алгоритм HI (или Н2), хранятся по строкам С, это будет специально оговариваться.
