Численное решение задач метода наименьших квадратов - Лоуcон Ч.
Скачать (прямая ссылка):
форму Ас =*у, где ац = qf(xt), / = 1.....т, / = 1,..., л + 2. Здесь с
л + 2-вектор с компонентами с^ а у - m-вектор с компонентами у у
Пусть заданные точки упорядочены так, что xt <... <хт. В S существуют базисы с тем свойством, что соответствующие им матрицы А ленточные с шириной ленты 4. Методы вычисления таких базисных функций для неравноудаленных узлов описаны в [36, 41, 45, 46]. Чтобы избежать усложнений, не существенных для иллюстрации ленточной последовательной
Таблица 27.1
Исходные данные для сглаживания посредством кубических сплайнов
2 4 6 8
10
12
2,2 4,0 5,0 4,6 2,8 2,7
14 16 18 20 22 24
3,8 5,1 6.1 6,3 5,0 2.0
Таблица 27.2
RMS как функция от NBP
nbp 5 6 7 В 9 10
RMS 0,254 0,085 0,134 0,091 0.007 0,0
172
обработки, мы ограничимся здесь обсуждением только случая равноудаленных узлов.
Чтобы построить нужный базис, положим h = bk+l — bk, k = l,...,n — I. Введем два кубических многочлена:
p,(r) = 0,25r3, p2(r)= 1 -0,75(1 +r)(l - г)2.
Пусть 11 обозначает замкнутый интервал [bit 03], а 1к — полуоткрытый интервал (bk,bk+l ],к = 2,... ,и - 1.
В интервале 1к только четыре из функций qf принимают ненулевые значения. Эти четыре функции определяются для х 61к формулами х - Ьк
Чк(х)=р1(1-t), <7*+i(jO=/>2(1 -t).
Рассмотрим иллюстративную задачу сглаживания, используя эту систему базисных функций и применяя к полученной ленточной матрице последо-
173
вателыюе хаусхолдерово накапливание. Исходные данные для этого примера приведены в табл. 27.1.
Число узлов указывает параметр NBP. Ему последовательно были присвоены значения 5, 6, 7, 8, 9, 10. При фиксированном значении NBP абсциссы узлов определяются формулами
22(/-1)
й, = 2 +-*-/-1.....NBP.
' NBP-1
Количество вычисляемых коэффициентов равно NC = NBP + 2. Заметим, что в случае NBP = 6 число коэффициентов совпадает с числом заданных точек. В этом случае сглаживающая кривая интерполирует эти. точки. Положим
RMS :
/1 12 Л1
где г i — невязка в /-й точке. Значение RMS для каждого случая указано в табл. 27.2. Отметим, что RMS не является монотонной функцией от NBP. График каждой из шести сглаживающих кривых приведен на рис. 27.6.
§ 5. Удаление строк
Опишем три метода для удаления строки из задачи наименьших квадратов Ах 3s Ь. Удобно будет ввести обозначение
С~[Л:Ъ]. (27.31)
Пусть С имеет т строк, и столбцов и ранг к. Предположим, что для С вычислен множитель Холесского R. Здесь R - верхняя треугольная А: X «-матрица такая, что
ос-[о] (27-32>
для некоторой ортогональной матрицы Q порядка т. Матрица/! удовлетворяет еще соотношению
RTR=CTC. (27.33)
Удобно сделать дополнительное предположение, что все к диагональных элементов R ненулевые. Так обязательно будет, если А: = и. Если же к < и, то невырожденности диагональных элементов можно добиться, переставляя в случае необходимости столбцы и проводя новую триангуляризации).
Пусть vT - строка С, которую нужно удалить. Без потери общности можно считать, что vT - последняя строка С. Пусть С — подматрица, образованная прочими строками С: [
C = [CJ. (27.34)
Заметим, что ранг С равен к либо к - 1. Практические приложения операции удаления данных чаще всего связаны только со случаем, когда
174
rank С = rank С = л. Мы рассмотрим, однако, общий случай rank С - 1 < rank С < rank С = к < и
в той мере, в какой он не потребует значительных дополнительных умозаключений.
Нужно найти множитель Холесского R для С, т.е. верхнюю треугольную к X п-матрицу того же ранга, что и С, для которой
ОС
где Q — некоторая ортогональная матрица порядка т летворяет еще соотношениям
RTR~=CTC = CTC- w
T=RTR- wT.
(27.35) 1. Матрица R удов-(27.36)
Первые два метода, которые мы опишем, были предложены в [68].
1. Удаление строки. В этом методе предполагается, что известна не только матрица R, но и матрица Q из (27.32). Представим Q в.блочном виде и перепишем формулу (27.32) :
м Г Qi Р ' [с 1 R
1{ иТ а • \ А - 0
т - к - 1 { Я . Lv J .0 .
(27.37)
тп
1
1
Умножим обе части (27.37) слева на ортогональную матрицу, которая преобразует вектор q в нулевой, изменяет а и не меняет первые к строк Q. В качестве такой матрицы можно взять, например, одно преобразование Хаусхолдера или произведение m - к - 1 преобразований Гивенса. Равенство (27.37) приобретает вид
б.
л т
и
-Л
0.г
Р
л
а
0j
(27.38)
Теперь умножим обе части (27.38) слева на последовательность к преобразований Гивенса, которые постепенно аннулируют элементы р, меняя при этом о. Пусть G'f обозначает матрицу Гивенса^ оперирующую со строками / и /'. Первое из левых умножений будет на Gkk+l, второе - на &к-1, *+i и тд. Последнее умножение будет на G, (к+1. Эта последовательность операций обеспечивает, что преобразованная матрица R, которую мы обозначим R, сохраняет верхнюю треугольную форму.
После всех преобразований равенство (27.38) переходит в
LQ2 0j
[дА-
W
О
(27.39)
Так как а - единственный ненулевой элемент своего столбца, а столбец
175
имеет единичную евклидову длину, то а = ± 1. Поскольку евклидовы длины строк также равны единице, то й должен быть нулевым вектором. Это
