Г.М. Кобельков - Курс лекций по численным методам (1160467), страница 3

Просмтор этого файла доступен только зарегистрированным пользователям. Но у нас супер быстрая регистрация: достаточно только электронной почты!

Текст из файла (страница 3)

Дискретное преобразование ФурьеСопоставим функции f : [0, 1] → R из пространства L1 её ряд Фурье:Xf (x) ∼ak e2πikx .(21)k∈ZБудем считать, что функция такова, чтоXk(22)|ak | < ∞.Зафиксируем некоторое число N ∈ N. Положим h := N1 . Рассмотрим равномерную сетку Ωn из N узлов наотрезке [0, 1], т. е.Ωn := {x = n · h, n = 0, . .

. , N − 1} .(23)Пусть k = mN + j, где 0 6 j < N , а x = nh ∈ Ωn . Тогда можно переписать ряд в видеXk∈Zak e2πikx=Xk∈Zak e2πik(nh)=−1X NX!amN +j e2πi(mN +j)nh =m∈Z j=0!=N−1XXj=0 m∈ZamN +j e2πi(mn) 2πij(nh)e=N−1 XXj=0m∈ZamN +j e2πijx . (24)В переходе, отмеченном «!», мы воспользовались абсолютной сходимостью, а потом — тем, что N h = 1 иe2πi(mn) = 1.ОбозначимXfj :=amN +j .(25)m∈Z8Таким образом, если ряд для функции f сходится абсолютно в точках сетки Ωn , то для x ∈ Ωn имеет местоформула дискретного преобразования Фурье:f (x) =N−1Xfj e2πijx .(26)j=0Ясно, что эта формула будет работать и для функций, заданных лишь на дискретном множестве Ωn : для этогонужно взять функцию, задать её в узлах, а далее продолжить константой на каждый отрезок и для неё ужевзять разложение Фурье.2.2.2.

Быстрое дискретное преобразование ФурьеНайдём коэффициенты fk разложения функции f по экспонентам. Рассмотрим дискретное скалярное произведениеN−1X(f, g) :=hf (nh)g(nh),(27)n=0где h — расстояние между узлами.Утверждение 2.3. Функции вида ek (x) := e2πikx ортогональны при разных k. Имеем X 2πi(k−l)(nh)(ek , el ) = e2πikx , e2πilx =he.(28)Если k = l, то (ek , el ) = 1, а если k 6= l, то получается геометрическая прогрессия, в формуле для которой вчислителе стоит разность двух экспонент с показателями, кратными 2πi. Значит, при k 6= l имеем (ek , el ) = 0. У нас есть формулы для прямого и обратного преобразования Фурье:fk = (f, ek ) =N−1Xhf (nh)e−2πik(nh) ,f (nh) =n=0N−1Xfk e2πik(nh) .(29)k=0Теперь мы будем решать такую задачу: по известным коэффициентам fk хотим определить значения функции в узлах, то есть числа f (nh).Если тупо вычислять значения коэффициентов, то мы получим сложность O(N 2 ).

Это очень медленно и напрактике не применимо. Сейчас мы получим гораздо более быстрый алгоритм, работающий за время O(N log N ).Пусть N = N1 · N2 . Разложим число k по модулю N1 : пусть k = k1 N1 + k2 , где 0 6 k2 < N1 . Разложим такжечисло n по модулю N2 : пусть n = n1 N2 + n2 , где 0 6 n2 < N2 . Тогда имеемf (nh) =N−1Xk=0fk e2πik(nh) =NX2 −1 N1 −1Xfk1 N1 +k2 e2πi(k1 N1 +k2 )(n1 N2 +n2 )h =k1 =0 k2 =0NX2 −1 N1 −1Xfk1 N1 +k2 e2πik1 N1 n2 h e2πik2 nh , (30)k1 =0 k2 =0так как произведение скобок в показателе экспоненты равноk1 n1 N1 N2 h +k1 n2 N1 h + k2 (n1 N2 + n2 )h.| {z }| {z }(31)n1Отсюда получаем, что (30) преобразовывается к видуf (nh) =NX1 −1k2 =0NX2 −1fk1 N1 +k2 e2πik1 N1 n2 hk1 =0!e2πik2 nh .(32)Заметим, что выражение в скобках не зависит от n1 , а зависит только от n2 .

Значит, его можно вычислитьодин раз и подставлять в формулу N1 раз при вычислении тех значений f (nh), для которых n (mod N2 ) = n2 .Итого получаем N · N2 операцийи N · N1 операцийна вычисление функции.√ на вычисление коэффициентов,√√Рассмотрим случай N1 = N2 = N . Тогда будет O(N · N ) операций.

Однако N — это еще не обещанныйlog N . Но никто не мешает применить к внутренней сумме ту же процедуру, и так далее, пока не останется однослагаемое. Пусть N = 2l . Тогда всего будет O(N · l) = O(N log N ) операций.92.3. Разделённые разности2.3.1. Определение разделённой разности и её простейшие свойстваОпределение.

Разделённой разностью функции f нулевого порядка назовём саму функцию f (x). Разделённой разностью порядка 1 называется функцияf (x1 ; x2 ) =f (x2 ) − f (x1 ).x2 − x1(33)Разделённая разность порядка n определяется индуктивно через разделённую разность порядка n − 1:f (x1 ; x2 ; . . . ; xn ) =f (x1 ; . . . ; xn−1 ) − f (x2 ; . . . ; xn ).x1 − xn(34)Теорема 2.4.

Имеет место формулаf (x1 ; . . . ; xn ) =nXj=1Qi6=jf (xj ).(xj − xi )(35)Докажем индукцией по порядку разделённой разности. Имеем1x1 − xnn−1nXXf (xj )f (xj )f (x1 )f (xn )−+ n−1+ средние члены.= Qnnn−1QQQj=2j=1(x−x)(x−x)(xj − xi )ji 1i(xn − xi )i=1i6=ji=2i6=ji=2(36)i=1Приводя средние члены к общему знаменателю, легко1 видеть, что получилось то, что надо. Следствие 2.1.

Разность многочлена Лагранжа и приближаемой функции выражается через разделённыеразности порядка n. Совсем легко видеть, чтоLn (x) =nXj=1f (xj )ωn (x)Q.(x − xj ) (xi − xj )(37)i6=jТогда, умножив и разделив разность f (x)−Ln (x) на ωn (x), а заодно поменяв знак у скобки (x−xj ) в знаменателе,получаемnXf (xj ) f (x)Q+f (x) − Ln (x) =  ωn (x) = f (x; x1 ; . . . ; xn )ωn (x),ωn (x) j=1 (xj − x) (xj − xi )(38)i6=jчто и требовалось. Вспоминая формулу для погрешности приближения многочленом Лагранжа, получаем, чтоf (x; x1 ; . .

. ; xn ) =f (n) (ξ).n!(39)2.3.2. Интерполяционная формула НьютонаПусть нам даны точки x1 , . . . , xm+1 и заданы значения функции f (x) в этих точках. Рассмотрим также многочлен Лагранжа по точкам x1 , . . . , xm (просто выкинем последний узел). Тогда, очевидно, многочленD(x) := Lm+1 (x) − Lm (x) имеет нули в точках x1 , . . . , xm , а кроме того, deg D = m.

Значит, он пропорционален многочлену ωn (x). Пусть D(x) = Am ωm (x). Теперь вычислим константу Am . Применим формулу (38) кмногочлену Lm+1 , приближая его многочленом Lm и подставляя x = xm+1 :Lm+1 (xm+1 ) − Lm (xm+1 ) = Am ωm (xm+1 ) = Lm+1 (x1 ; . . . ; xm+1 )ωm (xm+1 ).(40)Осталось заметить, что Lm+1 (x1 ; .

. . ; xm+1 ) = f (x1 ; . . . ; xm+1 ) (потому что разделённая разность зависит толькоот значений в узлах), поэтому в общем случае имеем Am = f (x1 ; . . . ; xm+1 ).1Вкнижке Богачёва всё это тупо и подробно расписано на двух страницах. Если интересно — посмотрите там.10Теперь представим многочлен Ln (x) в виде(41)Ln = L1 + (L2 − L1 ) + (L3 − L2 ) + . . . + (Ln − Ln−1 ).Применяя предыдущие рассуждения, получаем интерполяционную формулу Ньютона:Ln (x) = f (x1 ) + f (x1 ; x2 )ω1 (x) + . . . + f (x1 ; .

. . ; xn−1 )ωn−1 (x).(42)Из этой формулы видно, как можно ускорить вычисление многочлена Лагранжа. Мы видим, что нам простонужно вычислить набор разделённых разностей. Будем вычислять их по такой схеме:(x(x111111))))ffff(x(x(x(x111111;;;;xxx222222))))ffff(x(xxff(x(xff(x222222))))(x(x(x111111;;;;xxffff(x(xxx222222;;;;xxxx333333))))(x(x222222;;;;xxffff(x(xxx333333))))ff(x(xff(x333333))))(x(x(x111111;;;;xxffff(x(xxx222222;;;;xxxx333333;;;;xxxx444444))))(x(x222222;;;;xxffff(x(xxx333333;;;;xxxx444444))))(x(x333333;;;;xxffff(x(xxx444444))))ff(x(xff(x444444))))(x(x(x111111;;;;xxffff(x(xxx222222;;;;xxxx333333;;;;xxxx444444;;;;xxxx555555))))(x(x222222;;;;xxffff(x(xxx333333;;;;xxxx444444;;;;xxxx555555))))(x(x333333;;;;xxffff(x(xxx444444;;;;xxxx555555))))(x(x444444;;;;xxffff(x(xxx555555))))ff(x(xff(x55555))))(x(x(x222222;;;;xxffff(x(xxx333333;;;;xxxx444444;;;;xxxx555555;;;;xxxx666666))))(x(x333333;;;;xxffff(x(xxx444444;;;;xxxx555555;;;;xxxx666666))))(x(x44444;;;;xxffff(x(xxx55555;;;;xxxx66666))))(x(x555555;;;;xxffff(x(xxx666666))))ff(x(xff(x666666))))(x(x(x111111;;;;xxffff(x(xxx222222;;;;xxxx333333;;;;xxxx444444;;;;xxxx555555;;;;xxxx666666))))ff(x(xxx222222;;;;xxxx333333;;;;xxxx444444;;;;xxxx555555;;;;xxxx666666;;;;xxxx777777))))ff(x111111;;;;xx(x(x(x222222;;;;xxffff(x(xxx333333;;;;xxxx444444;;;;xxxx555555;;;;xxxx666666;;;;xxxx777777))))(x(x33333;;;;xxffff(x(xxx44444;;;;xxxx55555;;;;xxxx66666;;;;xxxx77777))))(x(x444444;;;;xxffff(x(xxx555555;;;;xxxx666666;;;;xxxx777777))))(x(x555555;;;;xxffff(x(xxx666666;;;;xxxx777777))))(xfff(xxxf(x(x666666;;;; xx777777))))ff(x(xff(x777777))))(xНа это уйдёт по 2 вычитания и по одному делению на каждую разность.

Итого 3 операции. А всего этих22разностей n2 + O(n). Значит, всего будет 3 · n2 = 32 n2 операций. Экономия грошовая, но всё равно приятно.2.3.3. Интерполяция с кратными узламиРанее узлы xi были однократными, то есть нам были заданы только значения функции в этих узлах. Нотеперь мы хотим большего: пусть у нас заданы значения производных в этих узлах до порядка mi − 1 включительно (говорят, что узел xi имеет кратность mi ). Сведём задачу к предыдущей: «расклеим» узлы, рассмотревkε, где i = 0, mi − 1. Для таких узлов построим обычный многочлен Лагранжа.

Мы знаем,точки xek := xi + miчто f (x; x1 ; . . . ; xn ) =f (n) (ξ)n! .Подставляя точки xek и устремляя ε к нулю, замечаем, что ξ → xi , т. е.f (ex0 ; . . . ; xemi −1 ) →f (mi ).mi !(43)Это даёт повод для введения следующего определения:Определение. Разделённой разностью кратности m1 называется величинаf (x1 ; . . . ; x1 ) :=| {z }m1f (m1 −1) (x1 ).m1 !(44)Определение разделённых разностей для случая многих узлов с кратностями аналогично.Формулу для многочлена Лагранжа с кратными узлами проще всего выписать, расклеив узлы, а потомсклеив их заново. При этом получится такое выражение:Lm1 +...+mn (x) = f (x1 ) + f (x1 ; x1 )(x − x1 ) + . .

. + f (x1 ; . . . ; x1 )(x − x1 )m1 −1 +| {z }m1(45)+ f (x1 ; . . . ; x1 ; x2 )(x − x1 )m1 + f (x1 ; . . . ; x1 ; x2 ; x2 )(x − x1 )m1 (x − x2 ) + . . .Её можно получить так: выписать обычную формулу (считая узлы однократными), а потом склеивать множители в степени.112.4. Наилучшее приближение в нормированных пространствах2.4.1. Общая теорияПусть R — линейное нормированное пространство (ЛНП).

Рассмотрим f ∈ R.PПусть g1 , . . . , gn ∈ R, причёмgi линейно независимы. Мы хотим приблизить вектор f линейной комбинациейci g i .PОпределение. ПоложимP ∆ = inf kf − cj gj k. Если существует вектор c = (c1 , . . . , cn ), на котором этот infдостигается, то элементcj gj называется элементом наилучшего приближения.PВведём обозначение Ff (c) = kf − cj gj k.Теорема 2.5. Элемент наилучшего приближения всегда существует. Покажем, что Ff непрерывна. В самом деле, X XXX1 2 f −cg−−cg(c1j − c2j )gj 6|cj | · max kgj k .fj jj j 6 Рассмотрим F0 (c) на сфере S n . Поскольку сфера — компакт, функция F0 достигает там своего минимума.Положим Fe = minF0 (c). Заметим, что Fe 6= 0, потому что иначе функции gi были бы линейно зависимыми.Snkc1Заметим, что F0 (c) = |c| · F0 |c|, .

. . , c|c|n > |~c| · Fe . Зафиксируем число γ > 2kfe . Положим B = {c : |c| 6 γ}.FПоскольку шар — это компакт, то на нём функция Ff достигает своего минимума. Положим F 0 = min Ff (c).BИмеем F 0 6 Ff (0, . . . , 0) = kf k. Теперь покажем, что вне этого шара значения функции Ff заведомо далеки отминимума. ИмеемX2 kf k · FeFf (c) > cj gj − kf k > |c| · Fe − kf k >− kf k = kf k > F 0 .FeЭто означает, что вне шара функция заведомо не достигает минимума. Значит, она достигает его внутри шара Bлибо на его границе. Замечание.

Характеристики

Список файлов лекций