1611688890-f641c9ec8276824e4686da772eb56520 (826652), страница 17
Текст из файла (страница 17)
. . , n − 1.Чтобы восстановить S по S ′′ , нужно теперь взять дважды первообразную (неопределённый интеграл) от S ′′ (x). Выполнив два раза интегрирование (2.50), получим для x ∈ [xi−1 , xi ]S(x) = γi−1(xi − x)3(x − xi−1 )3+ γi+ C1 x + C26hi6hi(2.51)с какими-то константами C1 и C2 . Но нам будет удобно представитьэто выражение в несколько другом виде:S(x) = γi−1(x − xi−1 )3(xi − x)3+ γi+ K1 (xi − x) + K2(x − xi−1 ), (2.52)6hi6hiгде K1 и K2 — также константы.13 Насколько законен переход к такойформе? Из сравнения (2.51) и (2.52) следует, что C1 и C2 должны бытьсвязаны с K1 и K2 посредством формулC1 = −K1 + K2 ,C2 = K1 xi − K2 xi−1 .У выписанной системы линейных уравнений относительно K1 и K2определитель равен xi−1 − xi = −hi , и он не зануляется.
Поэтому переход от C1 и C2 к K1 и K2 — это неособенная замена переменных.Следовательно, оба представления (2.51) и (2.52) совершенно равносильны друг другу.13 Строго говоря, константы C , C , K , K нужно было бы снабдить ещё допол1212нительным индексом i, показывающими их зависимость от подинтервала [xi−1 , xi ],к которому они относятся. Мы не делаем этого ради краткости изложения.972.6.
СплайныДля определения K1 и K2 воспользуемся интерполяционными условиями. Подставляя в выражение (2.52) значения x = xi−1 и используяусловия S(xi−1 ) = yi−1 , i = 1, 2, . . . , n, будем иметьγi−1(xi − xi−1 )3+ K1 (xi − xi−1 ) = yi−1 ,6hiт. е.γi−1h2i+ K1 hi = yi−1 ,6откудаK1 =yi−1γi−1 hi−.hi6Совершенно аналогичным образом, подставляя в (2.52) значение x = xiи используя условие S(xi ) = yi , найдёмK2 =yiγi h i−.hi6Выражение сплайна на подинтервале [xi−1 , xi ], i = 1, 2, . . . , n, выглядитпоэтому следующим образом:S(x) = yi−1x − xi−1xi − x+ yihihi(xi − x)3 − h2i (xi − x)(x − xi−1 )3 − h2i (x − xi−1 )+ γi−1+ γi.6hi6hi(2.53)Оно не содержит уже величин αi , βi и ϑi , которые фигурировали висходном представлении (2.48) для S(x), но неизвестными остались γ1 ,γ2 , . .
. , γn−1 (напомним, что γ0 дано по условию задачи).Чтобы завершить определение вида сплайна, т. е. найти γ1 , γ2 , . . . ,γn−1 , можно воспользоваться условием непрерывности первой производной S ′ (x) в узлах x1 , x2 , . . . , xn−1 :S ′ (xi − 0) = S ′ (xi + 0),i = 1, 2, . .
. , n − 1.(2.54)Продифференцировав по x формулу (2.53), получим для x ∈ [xi−1 , xi ]S ′ (x) =3(xi − x)2 − h2i3(x − xi−1 )2 − h2iyi − yi−1− γi−1+ γi. (2.55)hi6hi6hi982. Численные методы анализаСледовательно, с учётом того, что xi − xi−1 = hi ,S ′ (xi ) ==yi − yi−1h23(xi − xi−1 )2 − h2i+ γi−1 i + γihi6hi6hihihiyi − yi−1+ γi−1 + γi .hi63(2.56)С другой стороны, сдвигая все индексы в (2.55) на единицу вперёд,получим для подинтервала x ∈ [xi , xi+1 ] представлениеS ′ (x) =3(xi+1 − x)2 − h2i+13(x − xi )2 − h2i+1yi+1 − yi− γi+ γi+1.hi+16hi+16hi+1Следовательно, с учётом того, что xi+1 − xi = hi+1 ,S ′ (xi ) ==3(xi+1 − xi )2 − h2i+1h2yi+1 − yi− γi− γi+1 i+1hi+16hi+16hi+1yi+1 − yihi+1hi+1− γi− γi+1.hi+136(2.57)Приравнивание, согласно (2.54), производных (2.56) и (2.57), которые получены в узлах xi с соседних подинтервалов [xi−1 , xi ] и [xi , xi+1 ],приводит к соотношениямhihi + hi+1hi+1yi+1 − yi yi − yi−1−,γi +γi+1 = 6 γi−1 +36hi+1hi(2.58)i = 1, 2, .
. . , n − 1,γ0 и γn заданы.Это система линейных алгебраических уравнений относительно неизвестных переменных γ1 , γ2 , . . . , γn−1 , имеющая матрицу2(h1 + h2 )h2h22(h2 + h3 )h31..,.h32(h3 + h4 )6.........00hn−12(hn−1 + hn )992.6. Сплайныв которой ненулевыми являются лишь три диагонали — главная и соседние с ней — сверху и снизу.
Такие матрицы называются трёхдиагональными. (см. §3.8). Кроме того, эта матрица обладает свойствомдиагонального преобладания (см. §3.2е, стр. 232): стоящие на её главной диагонали элементы 31 (hi + hi+1 ) по модулю больше, чем суммамодулей внедиагональных элементов этой же строки. В силу признакаАдамара (он рассматривается и обосновывается в §3.2е) такие матрицынеособенны. Как следствие, система линейных уравнений (2.58) относительно γi , i = 1, 2, .
. . , n−1, однозначно разрешима при любой правойчасти, а искомый сплайн всегда существует и единствен. Для нахождения решения системы (2.58) с трёхдиагональной матрицей может бытьс успехом применён метод прогонки, описываемый ниже в §3.8.Интересен вопрос о погрешности интерполирования функций и ихпроизводных с помощью кубических сплайнов, и ответ на него даётследующаяТеорема 2.6.1 Пусть f (x) ∈ Cp [a, b], p = 1, 2, 3, 4, а S(x) — интерполяционный кубический сплайн с краевыми условиями (I), (II) или (III),построенный по значениям f (x) на сетке a = x0 < x1 < .
. . < xn = bиз интервала [a, b], с шагом hi = xi − xi−1 , i = 1, 2, . . . , n, причём узлыинтерполяции являются также узлами сплайна. Тогда для k = 0, 1, 2и k ≤ p справедливо соотношениеmax f (k) (x) − S (k) (x) = O(hp−k ),x∈[a,b]где h = max hi .iПри формулировке этого утверждения и далее в этой книге мыпользуемся символом O(·) — «O-большое», введённым Э. Ландау и широко используемым в современной математике и её приложениях.
Длядвух переменных величин u и v пишут, что u = O(v), если отношениеu/v есть величина ограниченная в рассматриваемом процессе. В формулировке Теоремы 2.6.1 и в других ситуациях, где идёт речь о шагесетки h, мы всюду имеем в виду h → 0. Удобство использования символа O(·) состоит в том, что, показывая качественный характер зависимости, он не требует явного выписывания констант, которые должныфигурировать в соответствующих отношениях.Обоснование Теоремы 2.6.1 разбивается на ряд частных случаев, соответствующих различным значениям гладкости p и порядка производ-1002.
Численные методы анализаной k. Их доказательства можно увидеть, к примеру, в [11, 14, 32]. Повышение гладкости p интерполируемой функции f (x) выше, чем p = 4,уже не оказывает влияния на погрешность интерполирования, так какинтерполяционный сплайн кубический, т. е. имеет степень 3.
С другойстороны, свои особенности имеет также случай p = 0, когда интерполируемая функция всего лишь непрерывна, и мы не приводим здесьполную формулировку соответствующих результатов о погрешности(её можно найти, к примеру, в книге [11]).Отметим, что, в отличие от алгебраических интерполянтов, последовательность интерполяционных кубических сплайнов на равномерной сетке узлов всегда сходится к интерполируемой непрерывной функции. Это относится, в частности, к функции |x| из примера С.Н.
Бернштейна и к функции Υ (x) = 1/(1 + x2 ) из примера Рунге, рассмотренным выше в §2.5. Важно также, что с повышением гладкости интерполируемой функции до определённого предела сходимость эта улучшается.С другой стороны, интерполирование сплайнами иллюстрирует также интересное явление насыщения численных методов, когда, начинаяс какого-то порядка, увеличение гладкости исходных данных задачиуже не приводит к увеличению точности результата. Соответствующиечисленные методы называют насыщаемыми.
Напротив, ненасыщаемыечисленные методы, там, где их удаётся построить и применить, даютвсё более точное решение при увеличении гладкости исходных даных[40]. Основной недостаток понятий насыщаемости / ненасыщаемости состоит в трудности практического определения гладкости данных, которые присутствуют в предъявленной к решению задаче.2.6вЭкстремальное свойство кубических сплайновСплайны S(x), удовлетворяющие на концах рассматриваемого интервала [a, b] дополнительным условиямS ′′ (a) = S ′′ (b) = 0,(2.59)называются естественными или натуральными сплайнами. Замечательное свойство естественных кубических сплайнов состоит в том, чтоони минимизируют функционалE(f ) =Zab2f ′′ (x) dx,1012.6. Сплайнывыражающий в первом приближении энергию упругой деформациигибкой стальной линейки, форма которой описывается функцией f (x)на интервале [a, b].
Краевые условия (2.59) соответствуют при этом линейке, свободно закреплённой на концах.Как известно, потенциальная энергия изгибания малого участкаупругого тела пропорциональна квадрату его кривизны (скорости изгибания в зависимости от длины дуги) в данной точке. Кривизна плоской кривой, задаваемой уравнением y = f (x), равна, как показываетсяв курсах дифференциальной геометрии,f ′′ (x)1 + (f ′ (x))23/2(см. подробности в [37, 61]). Поэтому упругая энергия однородной линейки, принимающей форму кривой y = f (x) на интервале [a, b], приусловии приблизительного постоянства f ′ (x), пропорциональнаZ b2f ′′ (x) dx.aТеорема 2.6.2 Если S(x) — естественный интерполяционный кубический сплайн, построенный по узлам a = x0 < x1 < .
. . < xn = b, аϕ(x) — любая другая дважды гладкая функция, принимающая в этихузлах те же значения, что и S(x), то E(ϕ) ≥ E(S), причём неравенство строго для ϕ 6= S.Доказательство этого факта не очень сложно и может быть найдено,к примеру, в [2, 11, 35].Будучи предоставленной самой себе, упругая линейка, закреплённая в узлах интерполирования, принимает форму, которая, как известно из физики, должна минимизировать энергию своей упругой деформации. Таким образом, эта форма очень близка к естественномукубическому сплайну. Сформулированное свойство называют экстремальным свойством естественных сплайнов,14 , и оно служит началомбольшого и важного направления в теории сплайнов.Несмотря на красивые свойства естественных кубических сплайнов,следует отметить, что в реальных задачах интерполяции для получения наилучших результатов следует всё-таки пользоваться дополнительной информацией о производных интерполируемой функции —14 Иногдатакже говорят о вариационном свойстве естественных сплайнов.1022.