Амосов А.А., Дубинский В.А., Копченова Н.В. Вычислительные методы для инженеров (1994) (1095853), страница 64
Текст из файла (страница 64)
Понятие о дискретном преобразовании Фурье и тригонометрической интерполяции 1. Дискретное преобразование Фурье. В прикладных исследованиях широко используются различные варианты преобразования Фурье~ функций непрерывного аргумента, а также представление функций в виде сходящихся тригонометрических рядов (рядов Фурье). Известно, например, что всякая непрерывно дифференцируемая периодическая с периодом 1 функция ~ может быть разложена в ряд Фурье; ФО ~(х) = Е акехр(2хзЛх). (11.83) Л=:0 Здесь ~ — мнимая единица.
Коэффициенты разложения вычисляются по формулам 1 а~ = / Дх) ехр( — 2хйх) ~Ь о (11.84) У-1 Дх~) = Е а~ехр(2хйх), О 4~'с Х й=о (11.85) ~ Жан Батист Жозеф Фурье (1768 — 1830) — французский математик, один нз основоположников математической физики.
339 Однако во многих случаях функция ~ бывает задана лишь в конечном числе точек х~ — 7/Ф, ~' = О, 1,, У вЂ” 1. В этом случае аналогом формулы (11.83) является разложение вида ! Заметим, что это разложение имеет место тогда и только тогда, когда тригонометрический многочлен У-1 Яд(х) = Е аьехр(2хйх) а=о (11.86) интерполирует функцию / по ее значениям в точках х~, О < у < У. Выше (см. пример 11.2) было доказано, что система функций щ(х) = = ехр(2тйх), О < 1 < У ортогональна на множестве точек х =,у/У, О < 1' < У, причем (у~, явь) = У. Следовательно, разложение (11.85) действительно имеет место, причем в силу равенства (11.19) коэффициенты аь определяются по формуле 1 К-1 аь = — Е /(ц) ехр(-2хйхД, О 4 1 < У.
У 1.0 (11.87) Операцию преобразования набора значений /(х0), /(х1), ..., г (щ-1) в набор коэффициентов а0, а1, ..., ау1 принято называть ирямым диснрегннмм иреобраэованием Фурье, в обратную операцию — обратннм дискретнмм преобразованием Фурье. Осуществление этих операпдй является важной составной частью многих алгоритмов. Для удобства изложения введем обозначение и = ехр(2ю/У) и перепишем формулы (11.85), (11.87) в следующем виде: Ф-! /(х~) = Е щит, О <1< У, й=о (11.88) 1 К-1 а~ = — Е /(х~) ьг Ц О < 1 < У. У 1=о (11.89) 2. Быстрое дискретное преобразование Фурье. Если вычисления проводить непосредственно по формулам (11.88) и (11.89), то на выполнение каждого из преобразований потребуется примерно Уэ арифметических операций. (Здесь под арифметической операцией понимается умножение двух комплексных чисел с последующим сложением.
Величины М считаются вычисленными заранее.) Однако в случае, когда число У не является простым, количество арифметических операций, требуемых для вычисления по формулам (11.88) и (11.89), можно существенно уменьшить. Поясним сказанное на примере вычислений по формулам (11.88). (Вычисления по формулам (11.89) производятся аналогично с заменой м на ы1.) Пусть У = У1Уэ, где 2 < У1, 2 < Уэ — целые числа. Представим индекс 1' в виде у' = 11У1 + ~, где О < н < Уъ О ~ га < У1. Положим Й = 340 = Ц1Кг + ао, где 0 4 а1 < У1, 0 < «о < Уг. Пользуясь тем, что 1г = %л~ = $Д~У+ йДоУг + йод и м = 1, имеем оЖ = ш ~ы . Заменяя в формуле (11.88) суммирование по индексу й операцией повторного суммирования по индексам 4> и Й~, получим г 1 г у( ) = Е Е „ ,„ ' г М = Е И~,л)~1о~, (11.90) аоо а1о 1 г о во=о где 43Фг У вЂ” 1 а(Ц»,1о) = Е а),у,), ы а~о 1г о (11.91) 3 а м е ч а н и е 1.
Часто разложение (11.85) записывают в эквивалентном виде: 341 ~ч Массив а содержит У чисел и для его вычисления требуется УЮ1 М арифметических операций. После того как найдены значения а (ао, ьг), на вычисления по формуле (11.90) требуется УХг операций. Таким образом, общее число арифметических операций равно У(Х1 + Фг). Заметим, что достигнута зкономия в числе операций, поскольку Х, + + Уг < Ф1Фг = У, как только Х > 4. Выигрыша удалось достичь благодаря тому, что оказалось возможным выделить группы слагаемых (11.91), которые используются для вычисления значений ~(гз) при различных г, но сами вычисляются лишь однажды.
Указанная выше идея развита в алгоритмах быстрого дискрепьного преобразования Фурье. В случае, когда У = У1 ° Фг ... Юп, (где 2 < Уз), с помощью быстрого дискретного преобразования Фурье можно выполнить дискретное преобразование Фурье за У (У~ + Уг + ... + У„,) арифметических операций. Особенно эффективным является этот алгоритм, когда число У является степенью числа 2 (т. е.
Ф = 2п). В этом случае вместо ЛР операций требуется выполнить лишь 2У!о8гУ операций. Например, для У = 1024 = 2'о этот алгоритм позволяет У ускорить вычисления в — = 1024/20 н 50 раз. 21о8гУ Широкое внедрение алгоритма быстрого дискретного преобразования Фурье в практические вычисления привело к подлинной революции во многих областях, связанных с обработкой числовой информации. Программы, реализующие различные варианты этого алгоритма, входят в стандартное математическое обеспечение ЭВМ и доступны массовому пользователю. /(х~) = Е ау.хр(2чйх~), -Ф/2< МЖ/2 что соответствует ''интерполяции тригонометрическим многочленом Ял(з) = Е а~ехр(2хЫт) -ж/2< ~ЗЧг (11.92) Здесь коэффициенты аь по-прежнему задаются.
формулой (11.87). 3 а м е ч а н и е 2. Хотя интерполяционные тригометрические многочлены (11.86), (11.92) и совпадают в точках х~, они принимают существенно разные значения в точках х, отличных от узловых. 3. Тригонометрическая интерполяция. Рассмотрим кратко задачу интерполяции функции /, заданной в точках 0 < з0 < х1 « ...
зу1 1 1 тригонометрическим многочленом (11.92). К ней приводит, например, типичная радиотехническая задача о тригонометрической интерполяции периодического сигнала. Не вдаваясь в довольно сложную проблему оценки погрешности тригонометрической интерполяции, отметим тем не менее, что для гладкой периодической с периодом 1 функции / есть основание рассчитывать на выполнение приближенного равенства /(х) м Я~(,с) для всех г б [О, Ц. Рассмотрим важный вопрос о чувствительности многочлена Яу к погрешностям в исходных данных. Пусть значения у,'. ю Дх;) интерполируемой функции задаются с погрешностями е; и известно, что ~ я,~ 1 1 Ь(у') для ~ = О, 1, ..., У вЂ” 1. Тогда вычисляемый по значениям у,* тригонометрический интерполяционный многочлен Я„' содержит пог- решность.
Для нее справедлива оценка сХ(,Я) = шах ~Ях) — Яд(х) ~ 4 ЛуЬ(у ), [О, 1] и аналогичная оценке (11.61) для алгебраических многочленов. Здесь Лу — постоянная, являющаяся аналогом константы Лебега Лу. Примечательно то, что в отличие от задачи интерполяции алгебраическими многочленами (см. ~ 11.10) оптимальным (т, е. дающим РМ минимальное значение Лу) является равномерное распределение узлов„ и которому отвечает значение Лу м — 1п [(Л + 1)/2). 342 Таким образом, при тригонометрической интерполяции выбор узлов ху — — 11'У (О 1 1 < Ж) является наиболее естественным с точки зрения как простоты вычисления коэффициентов многочлена (быстрое дискретное преобразование Фурье), так и минимизации влияния ошибок исхОДных Данных.
$11.13. Метод наименьших квадратов Задача наименьших квадратов возникает в самых различных областях науки и техники. Например, к ней приходят нри статистической обработке экспериментальных данных с помощью регрессионного анализа, В инженерной деятельности задача наименьших квадратов используется в таких областях, как оценивание параметров и фильтрация. 1. Линейная задача наименьших квадратов. Пусть функция у = = 1 (х) задана таблицей приближенных значений у,м~(х), а=0,1,...,а, (11.93) полученных с ошибками е; = у, — вь где у; = ~(ю;), Если значения у~ о о получены из эксперимента, то ошибки носят случайный характер и зачастую уровень погрешности (" шума" таблицы) бывает значительным (рис.
11.11, а). Рис. 11.11 Предположим, что для аппроксимации функции 1 используется линейная модель: (11.94) у = Ф„,(х) = аорто(х) + еще~(х) + ... + а„уъ~х). Здесь ~ро(х), ~р~(х), ...„~р„,(х) — заданные базисные функции~, ао аь " В данном параграфе рассматриваются функции, .принимающие только вещественные значения. 343 аа — параметры модели, являющиеся коэффициентами обобщенного многочлена фа(х). Как уже отмечалось выше, одной из наиболее прос- тых и часто используемых линейных моделей вида (11.94) (при у~(х) = = х") является полиномиальная модель у = Р.
(*) - =о + а1э+ — + а ". (11.95) В случае, когда уровень неопределенности исходных данных высок, неестественно требовать от модели (11.94) выполнения условий (11.7) совпадения значений обобщенного многочлена ф (х) в точках х, с заданными значениями уь т. е. использовать интерполяцию. Как нетрудно видеть (см. рис. 11.11, б), при интерполировании происходит повторение ошибок наблюдений, в то время как при обработке экспериментальных данных желательно, напротив, их сглаживание. Отказываясь от требования выполнения в точках х; точных равенств (11.7), следует все же стремиться к тому, чтобы в этих точках выполнялись соответствующие приближенные равенства ОИэ( Ь) + о1й( 0) + - + Юв(~0) я уо, о~в( 1)+ а~я(~) + - + М 1) у1, (11.96) оЮэ(~) + оЮ( ) + - + М ) я у .
Используя обозначения (11.9), запишем систему приближенных ра- венств (11.96) в матричном виде: (11.97) Рва у. б(ф~, у) = обобщенного многочлена ф„(х) = Е а у.(х) от заданных табличных 1=О значений у, (О < 1 < л). Заметим, что минимум среднеквадратичного уклонения достигается при тех же значениях ао, а1, ..., а~, что'и минимум функции 344 Из различных критериев, позволяющих выбрать параметры ао, аь ..., а,„модели (11.94) так, чтобы приближенные равенства (11.96) удовлетворялись наилучшим в некотором смысле образом, наиболее часто используется критерий наименьших квадратов. Согласно этому критерию параметры выбираются так, чтобы минимизировать средмеквадраижчяое уклокеиие и а в(а, у) — Е Е (аурЯхг) — у,)) = $Ра — у$~ г О,г=о причем б (Ф у)' = ( у). 1 Итак, линейная задача метода наименьших квадрапгов состоит в следующем.
Требуется найти (при фиксированном наборе функций ро, гог, ..., уи) обобщенный многочлен Фй(х), для которого среднеквад- ратичное уклонение принимает минимальное значение: б (Ф$, у) = ш1п б (Ф„, у). Фи — =О, 1=0,1, ..., т. дв дж (11.98) Вычисляя частные производные функции в и изменяя порядок сум- мирования, от равенства (11.98) переходим к системе линейных алгеб- раических уравнений и п Е ( Е о,.(х,) уг(х;)) а; = 1, у,уь(х;) (х = О, 1, ..., т), у=о г=о г=о (11.99) которая называется нормальнои системой .иетода наименьших «вадратов.
Как нетрудно видеть, нормальную систему можно записать в виде (11.100) Р"Ра = Ргу, или, используя матрицу Грама Г = Р'Р (см. ~ 11.2) и вводя вектор б = = Р'у, в виде (11.101) (11.101). Тогда для Га= б. Л е и и а 11.1. Пусть °: решение системы любого а' = а + Ьа имеет место равенство в (а' у) = в (а, у) + 1РЬа12. Искомый обобщенный многочлен Фв будем далее называть многоиленом наилучшего среднеквадратичного приближения. 2. Нормальная система.