Численное решение задач метода наименьших квадратов - Лоуcон Ч.
Скачать (прямая ссылка):
и если только участвующие здесь матрицы не обладают очень специальной структурой, то машинное представление матрицы СТС + е2ЕтЕ три достаточно малом е будет неотличимо от представления матрицы СТС (а вектора СTd + e2ETf- от вектора CTd).
(CTC + e2ETE)x=CTd+e2ETf,
114
Тем не менее Пауэлл н Рид [146] обнаружили экспериментально, что приемлемое решение задачи с сильно различающимися весами типа задачи (22.1) можно получить преобразованиями Хаусхолдера, если позаботиться о введении в алгоритм специальных перестановок строк. Отсылаем читателя к гл. 17, где основные теоретические результаты работы [146] представлены в виде трех теорем. В названной работе построению k-го преобразования Хаусхолдера предшествуют следующие перестановки. Во-первых, как и в алгоритме HFTI (см. 14.9), выполняется обычная перестановка столбцов, переводящая в положение k-то столбец с наибольшей евклидовой длиной отрезка от к-й строки до m-й. Затем просматриваются элементы столбца А: с номерами к,..., т. Если элемент с наибольшим модулем находится в строке /, то строки Ink переставляются.
Для задачи (22.1) это правило перестановок будет большей частью чересчур перестраховочным. Главное соображение, лежащее в основе анализа Пауэлла и Рида, таково: если (в обозначениях формул (10.5)-(10.11)) главный элемент vp значим для задачи, но незначителен по величине в сравнении с некоторыми элементами и,-, i >р, то может произойти существенная потеря численной информации. Под "значим для задачи" мы подразумеваем, что ее решение сильно изменилось бы, если vp заменить нулем. Но если I vp | достаточно мал сравнительно с каким-либо числом | и,-1, i>p, то при машинном вычислении s и ир по формулам (10.5) и (10.7) vp будет неотличим от точного нуля.
Предположим, что ненулевые элементы матриц С и ? задачи (22.1) имеют одинаковый порядок, так что если и есть сильный разброс в порядках ненулевых коэффициентов матрицы этой задачи, то он связан только с малым параметром е. Тогда отмеченная выше катастрофическая потеря численной информации могла бы, как правило, произойти лишь в том случае, если бы главный элемент vp (из (10.5)) был выбран из строки матрицы еЕ (скажем, vp = ее,,-), хотя некоторая строка С (скажем, 1-я) еще не выбиралась в качестве ведущей и содержит такой элемент с,,-, что | СцА > е | eti |. Этой ситуации можно избежать, если сохранять порядок строк, указанный в (22.1), так, чтобы все строки С были использованы в качестве ведущих прежде, чем будет использована какая-либо строка еЕ. Здесь мы применили символы е? и С для обозначения как матриц исходной задачи, так н их наследниц в алгоритме.
Перейдем теперь к анализу сходимости решения задачи (22.1) к решению задачи НКУ при е -+ 0. Рассмотрим вначале специальный случай задачи НКУ, когда С - диагональная матрица, а ? - единичная матрица. Этот случай легко разобрать, а, как будет показано в дальнейшем, анализ более общих случаев сводится к этому специальному. Положим s, 0 0. . . 0 -j
,*,>...>5И1 > 0, эн(22.2)
J т , X п
? = /„. ' "(22.3)
Для специального случая (22.2), (223) приводимые ниже две леммы показывают, что при е->0 решение задачи наименьших квадратов (22.1) сходится к решению задачи НКУ. • :
8» 115
Лемма 22.4. Пусть Си Е - матрицы, определенные в (22.2) и (223). Тогда решение задачи НКУ в обозначениях (20.2), (203) запишется так:
xt =
, /=1,.
f(, i =m, + 1.....и.
При этом длина невязки Ex - f равна
(22.5)
(22.6)
Лемма 22.7. Пусть Си Е - матрицы, определенные в (22.2) и (22.3). Единственное псевдорешение задачи (22.1) запишется так:
f ¦¦<л~){'Н^)Т-'-1....."¦•<2г8)
i =т, +1,... ,л.
Лемму 22.7 легко проверить, построив и решив систему нормальных уравнений (CTC + e2ErE)x -CTd + e2ETf. При наших предположениях относительно С и ? эта система имеет диагональную матрицу коэффициентов, откуда легко следует (22.8).
Чтобы удобно записать разность между х(е) и х, определим векторно-значную функцию й(е) соотношением
х,(е) -х{ = е2И{(е), / = 1, .... и.
Согласно (22.5) и (22.8), sl4fi-d,ls,)
й,(е) =
1 +(в/1,У
0, i =m, +1.....л,
(219)
(2110)
|А/(е)1 <
»/ /|--
Is,* - ¦
0, i =m, + 1,... ,л. Вектор невязки для задачи (22.1) имеет вид
Г-е'СЙ(е) |
(22.11)
(22.12)
116
Введем взвешенную евклидову норму
Ме = [(-)2 + ...+(•)* +(./«)»+... + (. /е>«] ¦/», (22.13)
тогда
I /- I.
| | -[ ^] *(е)| =«4 ВСА(е) 1,2 + 1/-^*- е2?л(е) I2. (22.14) Из (22.10) видим, что х(е)-+х при е-*0, а из (22.14) следует, что при
'fa-см
При практических вычислениях представляет интерес вопрос, насколько малым следует взять е , чтобы гарантировать неразличимость х(е) и х, если их компоненты представлены машинными словами с относительной точностью rj. Fx л и все компоненты вектора х ненулевые, то это условие обеспечивается требованием, чтобы
|еаЛ,(«)1<ч1*/1. ''=1.....и- А (22.15)
Оно будет удовлетворено при всех I е | < «о. если J j Ь ys
[17s2 14/1,1 1 - I jjw
з . f 4*7 14
€o = mm {-
' ll/i-<
: /,-</,/«,* 0. j { (22.16)
Согласно (22.9) и (22.11), для нормы разности между x(f) я х можно дать оценку t ^ e2l/-xl
