Если у вас возникли сложности с курсовой, контрольной, дипломной, рефератом, отчетом по практике, научно-исследовательской и любой другой работой - мы готовы помочь.
Предоплата всего
Подписываем
Пусть требуется решить уравнение f(х) = 0, х[а,.b], (13.1)
т. е. найти все корни xi* удовлетворяющие этому уравнению на отрезке [а,b].
Задача численного решения уравнения f(х) =0 сводится, во-первых, к отделению корней, во-вторых, к последующему уточнению корней.
Отделить корни уравнения f(x)=0, х [а,b], значит заключить каждый корень xi* в интервал ai < xi* < bi , на котором имеется ровно одно решение.
Если функция f(х) такова, что без особого труда можно построить ее график, этим следует воспользоваться, чтобы представить ситуацию с количеством и расположением нулей f(х), выделяя те промежутки оси абсцисс, где график у=f(x) пересекает Ох. Описанный прием называют графическим способом локализации (другими словами, отделения, изоляции) корней.
Убедиться в том, что на данном отрезке [а,b] имеется нуль непрерывной
функции f(х), можно аналитическим способом, в основе которого лежит известное утверждение математического анализа.
Теорема 13.1 (Больцано – Коши). Если непрерывная на отрезке [а,b] функция
f(x) на концах его имеет противоположные знаки, т.е. f(a) f(b) <0, (13.2)
то на интервале (а, b) она хотя бы один раз обращается в нуль.
Результат, сформулированный в виде теоремы 13.1, можно значительно
усилить, если требование непрерывности функции f(x) на [а,b] дополнить
требованием монотонности ее на этом отрезке.
Теорема 13.2. Непрерывная строго монотонная функция f(x) имеет и притом единственный нуль на отрезке [а,b] тогда и только тогда, когда на его концах она принимает значения разных знаков.
Последняя теорема позволяет не только принимать, но и отвергать те или иные промежутки из области определения D(f) данной функции на предмет дальнейшего поиска ее нулей, если известно о ее монотонном поведении на этих промежутках и определены знаки значений функции на их концах.
Реально установить монотонность на данном отрезке можно для
дифференцируемой функции, потребовав знакопостоянство ее производной на всем отрезке. Для таких функций основой решения задачи локализации корней уравнения (13.1) может служить следующая теорема.
Теорема 13.3. Пусть f С1 [а,b]. Тогда, если f’(x) не меняет знак на (а, b), то
условие (13.2) является достаточным для того, чтобы уравнение (13.1) имело и притом единственный корень на отрезке [а,b].
Метод дихотомии:
Наиболее употребительным частным случаем метода дихотомии является метод половинного деления, реализующий самый простой способ выбора пробной точки — деление промежутка существования корня пополам. Выполнить приближенное вычисление с точностью корня уравнения f(х) = 0, х[а,.b], методом половинного деления при условии, что f(x) непрерывна на [а,b] выполнено условие (13.2), можно, например, по следующей схеме:
Шаг 0. Задать концы отрезка а и b, функцию f, малое число > 0 (допустимую абсолютную погрешность корня или полудлину его промежутка неопределенности), малое число > О (допуск, связанный с реальной точностью вычисления значений данной функции); вычислить (или ввести) f(a) и f(b). Шаг 1. Вычислить с := 0.5 (а + b). Шаг 2. Если bа<2, положить с ( корень) и остановиться. Шаг 3. Вычислить f(c). Шаг 4. Если | f(c)| <, положить b:=с и остановиться. Шаг 5. Если f(a)f(c) < 0, положить b := с и вернуться к шагу 1; иначе положить а:=с, f(b):=f(c) и вернуться к шагу 1.
Используемый в методе половинного деления способ фиксирования пробной точки можно охарактеризовать как пассивный, ибо он осуществляется по заранее жестко заданному плану и никак не учитывает вычисляемые на каждом шаге значения функции. Логично предположить, что в семействе методов дихотомии можно достичь несколько лучших результатов, если отрезок [а,b] делить точкой с на части не пополам, а пропорционально величинам ординат f(a) и f(b) графика данной функции f(x). Это означает, что точку с есть смысл находить как абсциссу точки пересечения оси Ох с прямой, проходящей через точки A(a,f(a)) и B(b;f(by), иначе, с хордой АВ дуги А%В (рис).
Запишем уравнение прямой, проходящей через две данные точки А и В: . Отсюда, полагая у = 0 (уравнение оси Ох), х = с (обозначение искомой точки пересечения прямой АВ с осью Ох), находим (13.3)
Метод, получающийся в развитие метода дихотомии таким фиксированием пробной точки, называют методом хорд, методом пропорциональных частей, методом линейной интерполяции.
4.Метод простых итераций решения нелинейных уравнений и систем. Теорема сходимости.
Для построения метода итераций преобразуем уравнение f(х) = 0, х[а,.b](13.1) к виду х= (х), х[а,.b]. Это можно сделать так: х= x+hf(x), тогда (х) = x+hf(x) . Постоянный множитель h выбирается при этом из условия сходимости метода (установим его позже). Пусть известно начальное приближение x0. Тогда x1 = (х0), x2 = (х1), … xn = (хn1). Приведенный способ построения числовой последовательности {xk} реализуется в методе итераций: xk = (хk1). k=l,2,... .
Рассмотрим, как ведет себя погрешность решения на итерациях метода. Обозначим хk+1 = х*+ k+1, хk = х*+ k, где k, k+1 — погрешности приближенного решения на двух соседних итерациях. Подставим представленные таким образом xk и xk+1 в итерационное правило: х*+ k+1 = (х*+ k). Разложим функцию в ряд Тейлора в окрестности точки х*: х*+ k+1 = (х*)+ k ‘(х*)+. Пренебрегая остаточным членом , получим соотношение, связывающее погрешность решения метода на двух соседних итерациях: |k+1 | = | k | | ‘(х*)|.
Сделаем некоторые выводы на основании этого приближенного равенства. 1)Если | ‘(х*)| < 1, то можно ожидать, что |k+1| < |k| и последовательность {xk} будет сходиться к решению. 2)Если | ‘(х*)| > 1, то скорее всего |k+1| > |k| и метод будет расходиться, так как каждое последующее приближение будет отстоять от решения дальше, чем предыдущее. 3)При 0< ‘(х*)<1 сходимость будет монотонной. 4)При 1< ‘(х*)<0 сходимость является немонотонной.
Теорема сходимости. Пусть уравнение (13.1) приведено к виду х= (х) таким образом, что функция (х) дифференцируема и выполняется условие | ‘(х)| < q<1 для всех х[а,.b].Тогда последовательность {xk} х*, причем ошибка
| xk х*| Замечание. Из теоремы следует, что для сходимости метода итераций необходимо, чтобы | ‘(х)| < q<1 для всех х[а,b]. Значит, |1+ hf'(x)| < 1. Отсюда можно сделать вывод, что при f'(x) <0 должно быть при f'(x) >0 должно быть
Метод итераций для системы уравнений. Пусть необходимо найти решение системы уравнений f1(х1, х2, …, хn) = 0, f2(х1, х2, …, хn) = 0, … fn(х1, х2, …, хn) = 0 (14.1). Применение метода итераций требует приведения этой системы к виду: х1=1(х1, х2, …, хn), х2=2(х1, х2, …, хn), … хn=n(х1, х2, …, хn).
В общем случае это можно сделать так же, как и для одного уравнения: хi = хi +h fi(х1, х2, …,); i = 1,2, .., n. Отсюда i(х1, х2, …, хn) = хi +h fi(х1, х2, …,). Параметр h здесь также выбирается из условия сходимости метода. Метод итераций для системы уравнений приобретает вид. Рассмотрим, как ведет себя погрешность на итерациях метода. Обозначим через вектор погрешности. Очевидно, что x(k) = x*+ (k), где х* точное решение. Как и для одномерного случая можно получить, что (k+1) = A (k), где
Теорема 14.1. Пусть система уравнений (14.1) приведено к виду х= (х) таким образом, что вектор-функция (х) дифференцируема и и все собственные значения матрицы A по модулю меньше единицы, т.е. |i |<1, i = 1,2, .., n, тогда последовательность {xk} х*.
5. Аналог метода Зейделя
Правило Зейделя.
Возвратимся к записи метода итераций для системы уравнений f1(х1, х2, …, хn) = 0, f2(х1, х2, …, хn) = 0, … fn(х1, х2, …, хn) = 0 в координатном виде. Сходимость итерационного процесса несколько улучшается, если при расчете функции i использовать вычисленные к этому моменту времени компоненты вектора x, а остальные компоненты вектора х взять с предыдущего шага:
6.Метод Ньютона для одного уравнения.
Рассмотрим уравнение f(х) = 0, х[а,.b], х* корень. Пусть найдено k-е приближение х(k). Представим х(k+1)= х(k)+х. Разложим в ряд Тейлора в окрестности точки х(k) функцию f(x)::
f(x) = f(х(k)+х) = f(х(k)) + f'(х(k))x + 0(x2).
Полагая, что x(k+1) = x* решение, ограничившись двумя членами разложения, получим
f(х(k)) + f'(х(k))x = 0.
Отсюда x = f(х(k))/ f'(х(k)), и итерационное правило метода Ньютона примет вид
х(k+1) = х(k) f(х(k))/ f'(х(k)), k = 1,2, ... . (14.3)
Здесь х(0) — начальное приближение.
Геометрический смысл правила Ньютона весьма прост. В плоскости ху построим график функции у = f(x)
Рис. 14.1
Точное решение х* будет точкой пересечения этого графика с осью абсцисс. Рассмотрим на этом графике точку Mk(х(k),f(х(k))) проведем касательную, проходящую через точку Mk. Точку пересечения касательной с осью х обозначим х(k+1).
Геометрическая интерпретация правила Ньютона может быть сформулирована следующим образом: последующее приближение х(k+1) ищется как точка пересечения с осью абсцисс касательной к функции f(x) в предыдущей точке х(k).
Теорема 14.2. Пусть функция f(x) удовлетворяет условиям
Тогда если члены последовательности {xk}, определяемые методом Ньютона (14.3), при любом фиксированном k = 1,2, ... принадлежат отрезку [а,b] и эта последовательность сходится на [а,b] к корню х* уравнения (13.1), то справедливы неравенства:
(Первое из этих неравенств позволяет считать (14.3) методом второго порядка, а второе, являясь апостериорной оценкой погрешности, может служить в качестве критерия останова процесса вычислений).
7.Видоизменения метода Ньютона. Метод Ньютона с кусочно-постоянной матрицей.
Перечислим недостатки метода Ньютона, на устранение которых направлены различные его модификации:
• трудность задания начального приближения, от которого метод сходится;
• необходимость вычисления на каждой итерации матрицы Якоби, что может потребовать существенных вычислительных затрат;
• необходимость решения на каждой итерации системы линейных алгебраических уравнений;
• требование невырожденности матрицы Якоби.
Рассмотрим модификации метода Ньютона, которые в той или иной мере устраняют перечисленные недостатки.
Метод Ньютона с кусочно-постоянной матрицей.Чтобы уменьшить вычислительные затраты на итерации, матрица Якоби в этой модификации остается постоянной на протяжении нескольких шагов. Число шагов т, на которых J постоянна, задается в такой модификации в качестве параметра, либо момент перевычисления матрицы Якоби определяется условием
в котором (0;l), например = 0.1 (матрица Якоби лишь при нарушении этого условия вычисляется заново).
Эффективность метода достигается в этом случае не только путем сокращения числа расчетов матрицы Якоби, но главным образом за счет того, что на т итерациях метода приходится решать линейные системы с одной и той же матрицей.
8. Видоизменения метода Ньютона. Методы продолжения по параметру
Эти методы позволяют обеспечить сходимость метода Ньютона от выбранного начального приближения x(0). Сущность методов заключается в замене исходной задачи последовательностью задач, так что каждая последующая задача при этом незначительно отличается от предыдущей. Последовательность строится таким образом, что первая система имеет решение, а последняя система совпадает с исходной задачей. Поскольку системы отличаются незначительно, то решение предыдущей задачи окажется хорошим начальным приближением для последующей. Решая такую последовательность задач методом Ньютона, получим в итоге решение исходной системы.
Пусть при решении системы F(x) = 0, x Rn.
используется начальное приближение x(0). Заменим исходное уравнение F(x) = 0 уравнением с параметром
H(x, t) = 0, t [0, 1],
которое при t = 0 имеет решение x(0), а при t = 1 совпадает с решением исходной задачи, т. е.
H(x(0), 0) = 0, H(x*, 1) = F(x*)= 0.
В качестве H(x, t) можно выбрать функции
H(x, t)=(t1) F(x(0))+ F(x)
либо
H(x, t)=( 1 t) (x x(0))+tF(x).
Разобьем отрезок [0, l] точками t0=0, t1, …, tm =1 на m интервалов. Получим искомую последовательность систем: H(x, ti)=0, i = 0,1, ..., m.
9. Метод Ньютона для систем нелинейных уравнений.
Пусть задана система нелинейных уравнений
f1(х1, х2, …, хn) = 0,
f2(х1, х2, …, хn) = 0, (14.1)
…
fn(х1, х2, …, хn) = 0.
решение которой достигается в точке x*= пространства Rn. Обозначим F=(f1, f2, … fn,))Т. Тогда исходная система запишется в виде F(x) = 0, х Rn.
Предположим, что известно k-e приближение х(k+1) к х*. Построим правило Ньютона вычисления (k+1)-гo приближения в форме х(k+1)= х(k) +х(k). (14.4)
Разложим функцию F(x) = F(х(k) +х) в ряд Тейлора в окрестности точки х(k) и сохраним в разложении два члена: F(x) = F(х(k) +х)= .
Полагая, что решение системы достигается на текущей итерации, относительно поправки х(k) получим систему линейных алгебраических уравнений: . (14.5)
Тогда и итерационное правило Ньютона решения системы нелинейных алгебраических уравнений запишется как
Такой вид метода Ньютона неудобен на практике, потому что требует вычисления обратной матрицы, а эта операция достаточно трудоемка. На практике метод Ньютона реализуется в следующем виде:
1. Решается система линейных алгебраических уравнений (14.5) и вычисляется вектор поправки.
2. Вычисляется (k+1)-e приближение по формуле (14.4).
3. Пункты 1, 2 повторяются для k = 1,2, ... до тех пор, пока не будет достигнута требуемая точность.
Критерием завершения итерационного процесса служат условия
или в более общей форме где числа и определены заранее.
10. Сведение решения системы нелинейных уравнений к решению вариационных задач.
Покажем, как поиск корня системы нелинейных уравнений можно свести к задаче минимизации и наоборот, минимизацию к решению систем уравнений.
Предположим, что в области D пространства Еn необходимо решить систему уравнений
f1(х1, х2, …, хn) = 0, f2(х1, х2, …, хn) = 0, fn(х1, х2, …, хn) = 0.
Пусть известно, что в D существует решение системы. Определим целевую функцию Ф(х) формулой
(1)
В области D целевая функция Ф(x) 0. Минимальное значение Ф(x) достигает система на векторе , который обращает все fi(х) в нуль. Поэтому решение системы эквивалентно поиску минимума Ф, определенной формулой (1): найти
Если минимальное значение Ф(x) в D строго больше нуля, то система уравнений не имеет решений. Если все функции fi(х) непрерывно дифференцируемы в открытой области D, то решить задачу можно методами анализа.
11. Метод покоординатного спуска
Идея всех методов спуска состоит в том, чтобы исходя из начального приближения точки перейти в следующую точку так, чтобы значение уменьшилось: Ф(х(1)) < Ф(х(0)).
Пусть в области D задано начальное приближение
Рассматриваем функцию при фиксированных значениях , как функцию одной переменной x. Находим
Значение x1, доставляющее минимум, обозначим :
Далее при фиксированных значениях рассматриваем как функцию одной переменной x2. Находим
Значение х2, доставляющее минимум, обозначаем , получаем и т. д., после n шагов получаем
В результате одного шага покоординатного спуска происходит переход из начальной точки в точку Если при этом оказывается, что
то начальный элемент x(0) точка, доставляющая экстремум Ф(x). Если Ф(х(1)) < Ф(х(0)), то выполняется следующий шаг покоординатного спуска, в котором за начальную точку принимается x(1), получаем x(2) и т. д. Этот процесс продолжается до тех пор, пока не выполнится какое-либо условие окончания алгоритма, например, |Ф(х(k+1)) Ф(х(k))| <, где заданная точность. (15.2)
12.Метод градиентного спуска.
Напомним, что градиент функции Ф(х) определяется формулой
Вектор gгad Ф(x) ортогонален линиям уровня Ф(х) = с = const, его направление совпадает с направлением максимального роста Ф(х) в заданной точке. В точке минимума функции gгad Ф(x) = 0. Определим итерационный процесс: х(k+1) = х(k) - h gгad Ф(x(k)), (15.3)
где h > 0 - шаг спуска, x(0) - заданное начальное приближение.
Заметим, что если gгad Ф(x) ¹ 0, то для достаточно малых значений h справедливо неравенство
Ф(x(0))> Ф(x(1)) > … > Ф(x(k)) > … .
Формула (15.3) представляет собой метод градиентного спуска с постоянным шагом h спуска. Итерационный процесс (15.3), например (15.2), либо следующего: ||gгad Ф(x(k))||< e,
Для примера рассмотрим поиск минимума функции Ф(х)= методом градиентного спуска. Согласно (15.3), имеем
Пусть начальное приближение , шаг спуска = 0,1. Первые три шага спуска лают следующие приближения
Увеличивать h можно только до определенного предела. Слишком большой шаг может привести к увеличению Ф(х). В рассматриваемом примере h = 2 приводит к и значению Ф(х(1)) = 9>Ф(х(0)) = 1,25. Таким образом, самым важным моментом в методе градиентного спуска является выбор шага. Различные стратегии выбора шага h=h(k) определяют различные варианты градиентного спуска.
Метод наискорейшего градиентного спуска.
Из (15.3) можно определить значение Ф(х) в точке x(k+1):
Если x(k) определено, то значение целевой функции Ф в следующей точке x(k+1) оказывается функцией только шага спуска. Будем выбирать h= из условия, чтобы функция Ф за этот шаг максимально уменьшила свое значение: .
Выбор шага hk сводится к отысканию минимума функции одной переменной.
В примере п. 15.3.2 выбор шага в точке х(0) сводится к следующей задаче:
.откуда определяется шаг h наискорейшего спуска h0 =10/9.
Заметим, что методы градиентного спуска, вообще говоря, определяют локальный минимум целевой функции Ф. Это связано с зависимостью всего пути спуска от начального приближения х(0). Известно, что различные начальные приближения х(0) дают пути, приводящие к различным точкам локального минимума. Поэтому очень важно в задача минимизации использовать всю имеющуюся информацию о зависимости целевой функции Ф от вектора х для правильного выбора начального приближения.
13. Постановка задачи интерполирования и ее разрешимость.
Пусть на отрезке [a, b] заданы n+1 точки a=x0, x1, …, xn=b, которые называются узлами интерполяции, и значения некоторой функции f(x) в этих точках f(x0)=y0, f(x1)=y1, …, f(xn)=yn (16.1)
Рисунок 16.1. Геометрическая интерпретация задачи интерполирования
Требуется построить функцию φ(x) близкую к f(x) и принимающую в узлах интерполяции те же значения, что и f(x), т.е. такую, что φ(x0)=y0, φ(x1)=y1, …, φ(xn)=yn.(16.2)
Функция φ(x) называется интерполирующей или интерполяционной для функции f(x) на отрезке [a, b], если ее значения в узлах интерполяции совпадают со значениями функции f(x).
Функция φ(x)должна обладать "хорошими" свойствами, позволяющими легко производить над ней те или иные аналитические или вычислительные операции. Геометрически это означает, что нужно найти кривую y=φ(x), проходящую через заданную систему точек Mi (xi , yi) (i=0, 1, 2, …, n) (рис. 16.1).
В такой общей постановке задача может иметь бесчисленное множество решений или совсем не иметь решений в зависимости от того класса функций, в котором ищется интерполяционная функции. Однако эта задача становится однозначной, если вместо произвольной функции φ(x) искать полином Pn(x) степени не выше n, удовлетворяющий условиям (16.2), т.е. такой, что Pn(x0)=y0, Pn(x1)=y1, …, Pn(xn)=yn. В этом случае будем говорить о полиномиальной аппроксимации. Для вычислительной математики многочлены привлекательны тем, что они являются линейными функциями своих параметров (коэффициентов), и их вычисления сводятся к выполнению конечного числа простейших арифметических операций – сложения и умножения.
Полученную интерполяционную формулу y=φ(x) обычно используют для приближенного вычисления значений данной функции f(x) для значений аргумента x, отличных от узлов интерполирования. Такая операция называется интерполированием функции f(x). При этом различают интерполирование в узком смысле, когда x[x0, xn], и экстраполирования, когда x[x0, xn]. В дальнейшем под термином интерполирование будем понимать как первую, так и вторую операцию.
14. Интерполяционный многочлен в форме Лагранжа
Интерполяционная формула Лагранжа используется для произвольно заданных узлов интерполирования.
Пусть в точках x0, x1, …, xn таких, что a x0<…< xnb, известны значения функции y=f(x), т.е. на отрезке [a, b] задана табличная (сеточная) функция f(x):
|
x |
x0 |
x1 |
… |
xn |
|
y |
y0 |
y1 |
… |
yn |
(1)
Нужно построить полином
Ln(x) =a0 + a1x + a2x2 + … + an xn (2)
степени n, имеющий в заданных узлах x0, x1, …, xn те же значения, что и функция f(x), т.е. такой, что
Ln(xi)=yi (i=0, 1, 2,…, n).
Это значит нужно найти его n+ 1 коэффициент a0, a1, a2 , …, an . Для этого имеется как раз п + 1 условие (φ(x0)=y0, φ(x1)=y1, …, φ(xn)=yn). Таким образом, чтобы многочлен (2) был интерполяционным для функции (1), нужно, чтобы его коэффициенты удовлетворяли системе уравнений
,Будем строить многочлен n-ой степени Ln(x) в виде линейной комбинации многочленов n-й же степени li(x) (i=0, 1, 2,…, n), где i - номер многочлена. Для того чтобы этот многочлен был интерполяционным для функции f(x), достаточно зафиксировать в качестве коэффициентов ci этой линейной комбинации заданные в табл. (1) значения yi=f(xi), а от базисных многочленов li(x) потребовать выполнения условия
(3)
где ij – символ Кронекера. В этом случае для многочлена в каждом узле xj (j{0, 1,…,n}) справедливо Ln(xj)=l0(xj)y0+…+lj-1(xj)yj-1+lj(xj)yj+lj+1(xj)yj+1+…+ln(xj)yn =0+…+ 0+yj+0+…+0=yj
li(x)=Ai(x-x0)…(x-xi-1)(x-xi+1)…(x-xn),
а коэффициент Ai этого представления легко получается из содержащегося в (3) требования li(xi)=1. Подставляя в выражение li(x) значение x=xi и приравнивая результат единице, получаем:
Таким образом, базисные многочлены Лагранжа имеют вид:
,
а искомый интерполяционный многочлен Лагранжа:
15.Остаток интерполирования в форме Лагранжа.
Многочлен Лагранжа: (1)
Будем выяснять величину отклонения f(x) от Ln(x) в произвольной точке x[a, b], иначе, величину остаточного члена Rn(x)= f(x) Ln(x)
интерполяционной формулы Лагранжа (1) в предположении, что f(x) Сn+1[a, b], т.е. данная функция n + 1 раз непрерывно дифференцируема.
Обозначим
Определ-й через узлы x0, x1, …, xn многочлен (n + 1)-й степени. Тогда можно показать, что (2)
где некоторая точка из промежутка интерполяции (a, b).
Теперь можно пытаться отвечать на вопросы о погрешности приближ-го вычисл-я значения f(x) с помощью Ln(x) в к-л конкретной точке промежутка [a, b] о величине max погрешности, допускаемой при подмене функции f(x) многочленом Ln(x) на отрезке [a, b] о сходимости интерпол-го процесса, т.е. о том, имеет ли место (f(x), Ln(x)) 0, при n по метрике (у) того или иного определ-го на [a, b] функцион-го пространства.
Так, если известна величина '
то оценить абсолютную погрешность интерполяц. формулы (1) в любой точке [a, b] можно с помощью неравенства
Max погрешность интерпол-я на отрезке [a, b] оценивается величиной
(3)
Так как максимумы модулей функций f (n+1)(x) и Пn+1(х) достигаются, вообще говоря, в разных точках отрезка [a, b], то более точной, но более трудно реализуемой по сравнению с (3) следует считать оценку
Недостатком интерпол-й формулы Лагранжа является то, что каждое слагаемое зависит от всех узлов интерполяции. При добавлении узла интерполяции и, следовательно, повышении порядка полинома необходимо вычислять не только слагаемое, относящееся к новому узлу, но и перевычислять заново слагаемые, относящиеся ко всем узлам интерполяции.
16.Разделенные разности и их свойства.
Понятие конечных и разделенных разностей
Пусть на равномерной сетке в узлах x0, x1, …, xn заданы значения функции f(x), т. е.
yi = y(xi), i = 0,…, n; хi = хi+1 + h; i == 0,…, n. h - шаг сетки.
Конечными разностями первого порядка называют величины: yi = yi+1yi, i = 0,…, n1.
Разности конечных разностей первого порядка называют конечными разностями второго порядка: Δ2yi= yi+1 yi, i = 0,…, n2.
Аналогично определяют конечные разности следующих порядков:
kyi = k1yi+1 k1yi, i = 0,…, nk.
Схема вычисления конечных разностей:
Если узлы интерполяции произвольные точки, то вводится понятие разделенных разностей.
Разделенная разность первого порядка:
Разделенная разность второго порядка:
Разделенная разность k-го порядка выражается через разделенные разности (k1) —го порядка, следующим образом
Свойства разделенных разностей:
1. Разделенная разность является симметричной функцией своих аргументов, т. е. значение разделенной разности не зависит от порядка перечисления узлов интерполяции, а зависит только от их значений.
2. Разделенная разность вычисляется и тогда, когда две или более точки интерполяции одинаковы. Разделенная разность первого порядка в случае, когда одинаковы две точки, имеет вид:
Если точка хi, среди узлов интерполяции встречается (m+1) раз, то разделенная разность m-ro порядка на этих узлах: 3. Если функция является многочленом степени n (), то ее разделенная разность n-го порядка постоянна, а разделенная разность (n+1)-го порядка обращаются в нуль.
Схема вычисления разделенных разностей:
17. Интерполяционная формула Ньютона для неравномерной сетки
Пусть х [a, b]. Построим первую разделенную разность
f(x,x0)=[f(x)-f(x0)]/(x-x0)
Откуда находим f(x)=f(x0)+(x-x0)f(x,x0)
Построим вторую разделенную разность: f(x,x0,x1)=[f(x,x0)-f(x0,x1)]/(x-x1)
И выразим из нее f(x, x0): f(x,x0)=f(x0,x1)+(x-x1)f(x,x0,x1).
Подставим f(x, x0) в выражение для f(x):
f(x)=f(x0)+(x-x0)f(x0,x1)+(x-x0)(x-x1)f(x,x0,x1)
Аналогично привлекаем следующие узлы интерполяции для построения интерполяционной функции. После использования всех узлов интерполяции:
F(x)=f(x0)+(x-x0)f(x0,x1)+(x-x0)(x-x1)f(x0,x1,x2)+…+(x-x0)…(x-xn-1)f(x0,…,xn)+(x-x0)…(x-xn)f(x,x0,…,xn)=Pn(x)+Rn(x).
Таким образом, интерполяционный многочлен Ньютона имеет вид:
Непосредственной проверкой нетрудно убедиться, что f(xk)=yk=Pn(xk), k=0,1,…n.
Погрешность интерполяции Rn(x)=(x-x0)…(x-xn)f(x,x0,xn).
Более эффективное вычисление значения функции по интерполяционной формуле Ньютона можно получить, если преобразовать ее к такому виду:
F(x)=f(x0)+(x-x0)(f(x0,x1)+(x-x1)(f(x0,x1,x2)+…+(x-xn-2)(f(x0,...,xn-1)+(x-xn-1)f(x0,…,xn)))…).
Интерполяционная формула Ньютона позволяет легко наращивать число узлов интерполяции, требуя при этом вычисления лишь дополнительных слагаемых. Например, добавления узла xn+1 приведет к вычислению слагаемого
(x-x0)…(x-xn)f(x0,x1,…,xn).
18.Конечные разности и их свойства. Интерполяционные формулы Ньютона для равномерной сетки. (КР=конечные разности)
Пусть на равном. в узлах x0, x1, …, xn заданы значения функции f(x), т. е.
yi = f(xi), i = 0,…, n; хi = хi+1 + h; i = 0,…, n.
КР 1-го порядка называют величины yi = yi+1yi , i = 0,…, n1.
Разности КР первого порядка назыв. КР 2-го порядка: 2yi = yi+1 yi, i = 0,…, n2.
Анал. опред. КР след.порядков kyi = k1yi+1 k1yi, i = 0,…, nk.Схему вычисления
x0 y0 y0
x1 y1 y1 2y0 3y0
x2 y2 y2 2y1 3y1 4y0
x3 y3 y3 2y2
x4 y4
Простейшие свойства КР:
1.КР постоянной =0;
2.постоянный мно-ль у ф-ции можно вынести за знак КР;ж
3.КР от суммы 2-х ф-ций =сумме их КР в одной и той же точке.
Интерполяционные формулы Ньютона.
Пусть ф-ция у = f(x) задана на сетке равноотстоящих узлов xi=x0+ih, где i = 0,1, ..., п, и для нее построена таблица КР.В соотв. с тем, что было сказано о направлении модификации интерполяц.ф-лы Лагранжа в начале предыдущего параграфа,будем строить интерполяц. многочлен Рп(х) в форме Рn(х) = а0+а1(хх0) + а2(хх0)(хх1)+... + аn(хх0)( хх1) … (ххn1). (1)
Его п+1 коэффициент а0, а1, ..., аn будем находить последовательно из п+1 интерпол. равенств Рn(хi) =yi, i = 0,1, ..., п.
А именно, полагая i = 0, т.е. х = х0, имеем Рn(х0) = а0, следовательно, а0= у0.
Далее, при i = 1 аналогично получаем равенство а0+а1(х х0)y1,
в кот. подставляем уже найденное значение а0= у0. Разрешая это ра-во относит. а1 и используя обозначение КР, получаем
Полной индукцией можно показать справедливость
Подставляя найденные коэф. а0, а1, ..., аn в (1), получаем многочлен
(2)
кот.назыв. первым интерполяционным многочленом Ньютона.
Учитывая, что каждое слагаемое многочлена (2), начиная со второго, содержит множитель х х0, естественно предположить, что этот многочлен наиболее приспособлен для интерполирования в окрестности узла х0.Узел х0 базовым для многочлена (2) и упростим (2) введением новой переменной q или равенством x = x0 +qh. Так как x xi = x0 + qh x0 ih= h (q i),то в результате подстановки этих разностей в (2) приходим к первой интерпол.формуле Ньютона в виде
(3)
где Рn(x0 + qh) указывает не только на n-ю степень многочлена, но и на базовый узел x0 и связь переменных х и q.
Первая формула Ньютона (3) применяется при |q|<1, а именно, для интерп-ния вперед ( х (х0, x1), т.е. при q (0, 1)) и экстраполирования назад (при х < х0 т.е. при q < 0).
В отличие от (1), форма интерполяц. многочлена Рn(х) берется такой, кот.предусматр.поочередное подключение узлов в обр. порядке: сначала последний, потом предпосл.и т.д.; .Р(х) = а0+а1(ххn) + а2(ххn)(ххn1)+... + аn(ххn)( ххn1)…(хх1).
Коэф.а0, а1, ..., аn этого многочлена нах-ся аналог.тому, как они находились для многочлена (1), только здесь подстановка узловых точек вместо х и рассмотрение интерпол. равенств производится тоже в обратном порядке.
Т.о., получаем второй интерполяционный многочлен Ньютона
(4)
в кот.базовым является узел хn и коэффициенты которого опр-ся КР, располож. на восходящей от уn диагонали.Положим в (4) xxn+qh, иначе, введем новую переменную и преобразуем к ней входящие в (4) разности:
x xi = xn + qh x0 ih= x0+ nh + qh x0 ih= h (q+n i)
В результате приходим ко второй интерполяционной формуле Ньютона вида
(5)
Ее целесообразно использовать при значениях |q| < 1, т.е. в окрестности узла хп для интерполирования назад (при q (1, 0)) и экстраполирования вперед (при q > 0).
Наряду со спец. вывед.для начала и конца таблицы 1-ой и 2-ой интерполяц. формулами Ньютона имеется еще неск.формул, рассчитанных на их применение в центральной части таблицы и потому назыв. центральными интерпол. формулами. Будем считать, что узел x0 расположен в середине таблицы, и нумерация ост.узлов производится, начинаясь с х0, с использованием как «+», так и «-«индексов, т.е. считаем xi=x0+ih, где i = 0, ±1, ±2,... . Все подчеркнутые в табл. КР наз-ся центральными разностями.
x3 y 3 y3
x2 y3 y2 2y3 3y3
x1 y1 y-1 2y2 3y2 4y3 5y3
x0 y0 y0 2y1 3y1 4y3 5y2 6y3
x1 y1 y1 2y0 3y0 4y1
x2 y2 y2 2y1
x3 y3
Для 1-го интерпол. многочлена Ньютона в форме (3) погрешность может быть записана след. образом
Для 2-го интерполяц. многочлена Ньютона в форме (5) погрешность может быть записана след. образом
19. Интерполяционная формула Стирлинга
Интерполяционный многочлен дан в форме
Р(х) = а0+а1(хх0) + а2(хх0)(хх1)+ а3( хх1) (хх0)(хх1)+
+а4( хх1) (хх0)(хх1)( хх2)+… .
Введя новую переменную и выразив через нее разности x xi = h (q i) для всех i= 0, ±1, ±2, ..., в результате подстановки этих разностей и выражений коэффициентов после преобразований приводит к формуле
называемой интерполяционной формулой Стирлинга.
Рассмотрим вопрос о том, как могут быть трансформированы остаточный член и его оценки при конечноразностной интерполяции.
Построенные здесь конечноразностный интерполяционные многочлены Стирлинга это всего лишь формы представления интерполяционного многочлена Лагранжа. Следовательно, для всех этих форм справедливо выражение остаточного члена:.
20. Многочлены Чебышева.
Определение. Многочленом Чебышева называется функция Тn (х):= cos(n arccos х), (1)
где n ∈ N0, x ∈ [-1,1].
Замечание. Многочлены Чебышева являются частными решениями дифференциального уравнения Чебышева при .
Убедимся, что функция Тn(х), представленная с помощью тригонометрических функций, на самом деле является многочленом при любом n∈N0.
Непосредственной подстановкой в (1) значений n=0 и n=1, получаем Т0 (х) =1,Т1 (х) = х.
Положив а:=arccos x, имеем: Т1 (х) =cos а, Тn (х)= cos nа, Тn−1 (х) = cos (n − 1)а, Тn+1 (х) = соs(n+ 1)а, и так как (по формуле суммы косинусов) соs(n + 1)а+cos(n − 1)а=2 cos a cos nа, то, значит, справедливо равенство Тn+1 (х) + Тn−1 (х) =2 Тn (х) Т1 (х). (2)
Формула (2) определяет при n =1, 2, 3,... последовательность функций Тn(х), начинающуюся с Т0 (х)=1, Т1(х)= х, рекуррентно; при этом нужно иметь в виду, что здесь
x ∈ [-1,1], как и в (1).
Подставляя в (2) заданные начальные члены последовательности {Тn(х)}, найдем несколько ее последующих членов:и т.д.
Анализ рекуррентной формулы (2) позволяет считать очевидными следующие факты:
1) все функции Тn(х), определенные в (1), являются многочленами при любом натуральном n;
2) степени этих многочленов возрастают с увеличением n, причем старший член многочлена Тn(х) равен
3) многочлены Тn(х) при четных n выражаются через степенные функции только четных степеней, при нечетных — только нечетных.
Наряду с многочленами Чебышева Тn(х) часто используют многочлены, получаемые из Тn(х) делением на старший коэффициент, т.е. — многочлены со старшим коэффициентом 1. Будем называть их нормированными многочленами Чебышева.
Многочлены Чебышева обладают рядом замечательных свойств. Рассмотрим некоторые их свойства, имеющие отношение к поставленной выше проблеме аппроксимации функций.
Свойство 1. Многочлен Чебышева Тn(х) имеет на отрезке [-1,1] ровно n различных действительных корней; все они задаются формулой где
Свойство 2. Корни многочленов Чебышева перемежаются с точками их наибольших и наименьших значений, равных соответственно +1 и −1 для Тn(х), а именно функция Тn(х), имеет экстремумы Тn(хj)= в точках
Свойство 3 (теорема Чебышева). Из всех многочленов степени n со старшим коэффициентом 1 нормированный многочлен Чебышева наименее уклоняется от нуля на отрезке [−1,1].
21 Минимизация остатка интерполирования
Если функция имеет непрерывную (n+1)-ю производную, то имеет место следующая оценка остаточного члена: , (1)
где , .
Погрешность интерполирования можно представить также через разделенную разность следующим образом:
Из формулы (1) следует, что для данной функции погрешность интерполирования зависит от выбора узлов на отрезке [a,b]. Величину можно минимизировать за счет выбора узлов интерполирования
, i=0,1,…,n.
Такими оптимальными узлами для отрезка [1,-1] являются корни многочлена Чебышева первого рода:
,
которые вычисляются по формуле:
, k=0,1,…,n.
В случае произвольного отрезка [a,b] из этого равенства получим формулу для оптимальных узлов:
, k=0,1,…,n.
При этом оценка (1.3) примет вид:
, и .
22. Интерполирование с кратными узлами
Пусть на промежутке [a, b] расположены m+1 несовпадающих узлов х0, х1,.., хm, и пусть в этих точках известны значения у0=f(х0), у1=f(х1),..., ym=f(xm) данной функции, а также некоторые ее производные(максимальный порядок производных в разных узлах различен; в каких-то узлах производные могут быть вовсе неизвестны). Такие узлы будем называть кратными узлами. Конкретнее, будем считать, что заданы:
в узле значения , …
в узле значения , …
…………………………………………………..
в узле значения , …
тогда кратность узла считается равной , узла — …
Предполагая, что суммарная кратность узлов есть + +... + =n + 1 ставим задачу построения многочлена степени n (не выше n), такого, что где m>0, значения функции f(x) и ее производных, и по определению считается . будем называть интерполяционным многочленом Эрмшпа, а совокупность условий—условиями эрмитовой интерполяции.
Формально можно считать, что нахождение такого многочлена состоит в том, чтобы однозначно определить n+ 1 коэффициентов а0, а1..., аn его канонического представления =
Совокупность вышеуказанных требований можно рассматривать как систему из п+1 уравнения относительно п+1 неизвестного — коэффициентов многочлена. Можно доказать существование и единственность многочлена , удовлетворяющего условиям эрмитовой интерполяции. Выявление общего вида интерполяционных многочленов Эрмита представляет непростую задачу и требует привлечения определенных сведений из теории функций комплексной переменной. Рассмотрим одну из возможных процедур фактического построения таких многочленов, не требующую знания их общего вида.
Пусть Lm(x) — интерполяционный многочлен Лагранжа, построенный по данным m+1 значениям f(хk)= yk, k= 0,1,..., m. Будем пользоваться обозначением Так как по условию число m заведомо не превосходит n, то по теореме о делении многочлена с остатком искомый многочлен Эрмита можно представить в виде =+ где — некоторый неизвестный пока многочлен степени n− m−1. Для построения многочлена будем привлекать информацию о производных данной функции, т.е. равенства в тех узлах хi, где первые производные заданы.
Продифференцировав равенство, имеем
Так как, то в тех узлах , где по условию эрмитовой интерполяции cправедливо можно записать Отсюда выражаем значения многочлена в узлах: = Правая часть этого равенства может быть вычислена; обозначим ее через . Таким образом, в ряде узлов xi, известны значения многочлена = , по которым этот многочлен однозначно восстанавливается обычной лагранжевой интерполяцией, если в условиях содержится только производные порядка выше первого. Если же в исходной информации имеются значения производных более высокого порядка, чем первый, то для восстановления многочлена ставится задача эрмитовой интерполяции, для чего, наряду с полученными его значениями , находят значения его производных путем дифференцирования равенства (возможно неоднократного, в зависимости от максимального порядка заданных производных функции). Эта процедура построения интерполяционных многочленов Эрмита все более низких степеней продолжается до исчерпания всей информации о функции и ее производных.
23-24. Понятие сплайн-функции. Сплайн-интерполирование. Построение кубического сплайна. Вариационная и физическая интерпретация кубического сплайна. Пусть на отрезке [а, b] задана упорядоченная система несовпадающих точек xk (k =0, 1, ..., n).
Опр. 1. Сплайном Sm(x) называется определенная на [а, b] функция, принадлежащая классу Сl[а, b] l раз непрерывно дифференцируемых функций, такая, что на каждом промежутке [xk1, xk] (k = l, 2,..., n) это многочлен n-й степени. Разность d = m l между степенью сплайна т и показателем его гладкости l называется дефектом сплайна.
Если сплайн Sm(x) строится по некоторой функции f(x) так, чтобы выполнялись условия Sm(xi) = f(xi), то такой сплайн называется интерполяционным сплайном для функции f(x). Примеры: 1. y = akx + bk с параметрами ak, bk, очевидно, является интерполяционным сплайном степени 1 дефекта 1.
Совпадение дефекта сплайна с его степенью обеспечивает просто непрерывность сплайна. Рассмотрим наиболее известный интерполяционный сплайн степени 3 дефекта 1. При этом будем исходить из предположения, что узлы сплайна а = х0, х1 х2, ..., хn= b одновременно служат узлами интерполяции, т.е. в них известны значения функции fk := f(xk), k = 0,1,..., n.
Опр. 2. Куб. сплайном дефекта 1, интерполирующим на отрезке [а, b] данную функцию f(x), называется функция
удовлетворяющая совокупности условий: а) g(xk) = fk (условие интерполяции в узлах сплайна); б) g(x)C2 [а, b] (двойная непрерывная дифференцируемость); в) g"(a) = g"(b) = 0 (краевые условия).
Физичечкая интерпретация: Определенный таким образом сплайн называют еще естественным или чертежным сплайном. Чертежники фиксировали в этих точках гибкую упругую рейку для построения непрерывной линии, которая соответствовала искомой функции, удовлетворяющей условиям а)-в).
(Вариационная): Достаточно показать, что функция g(x) доставляет минимум функционала, характеризующее величину потенциальной энергии закрепленной в (n + 1)-й точке упругой рейки. Для построения по данной функции f(x) интерполирующего ее сплайна нужно найти 4n его коэффициентов ak, bk, ck, dk (k = 1, 2, ..., n).
Имеем: из условий интерполяции а): g1(x0) = f0, gk(xk) = fk при k = 1, 2, ..., n;
из условий гладкой стыковки звеньев сплайна б): gk1(xk1) = gk(xk1), gk1(xk1) = gk(xk1),
gk1(xk1) = gk(xk1), при k = 2, 3,..., n; из краевых условий в): g1(x0) = 0, gk(xk) =0.
Как видим, условий оказалось 4n — ровно столько, сколько в записи сплайна (17.7) неизвестных коэффициентов. Подставляя сюда выражения функций и их производных через коэффициенты ak, bk, ck, dk при указанных значениях k и полагая для краткости hk= xk xk1 получаем детализированную систему связей
(17.8)
Теперь ставим задачу выявления эффективного способа нахождения коэффициентов сплайна ak, bk, ck, dk (k = 1,2, ...,n) из этой линейной относительно них системы (17.8). С этой целью будем исключать из системы неизвестные ak, dk, bk и сводить все к решению системы относительно неизвестных ck. Если ввести обозначение разделенной разности , то справедлива следующая
Теорема 17.1. При заданных в точках а и b краевых условиях (g (a) = g"(b) = 0) по заданным в узлах (17.6) значениям fk ( k = 1, 2,..., n,) функции f(x) можно построить единственный интерполирующий ее кубический сплайн g(x) дефекта 1, коэффициенты которого определяются системой
(17.13)
где k = 2,3,..., п, с0=0 и сn=0. Отсюда,
ak=fk.. Диагональное преобладание в трехдиагональной матрице системы (17.13) обеспечивает корректность и устойчивость метода прогонки:
при k = 3, 4,...,n, а затем к получению искомых значений ck обратной прогонкой по формуле
Теорема 17.2. Пусть g(x) — кубический сплайн дефекта 1, интерполирующий на системе узлов (17.6) отрезка [а, b] четырежды непрерывно диф-мую на нем функцию f(x). Тогда при любом фиксированном n найдется такая постоянная С > 0, что для любого х [а, b] справедливо неравенство:
25. Задача о наилучшем приближении в линейных нормированных пространствах.
Постановка задачи:
Задана функция y(x), некоторый класс функций X и некоторая норма, или метрика.
Найти функцию , такую что .
Чаще всего используются нормы и , такие что
(равномерное приближение);
(среднеквадратичное приближение).
Равномерное приближение:
Можно хотя бы частично ответить на вопрос о сходимости последовательности интерполяционных многочленов Лагранжа Ln(x) к функции f(x) на отрезке [а, b]. Для которой они построены: если функция f(x) бесконечно дифференцируема на [а, b] и в качестве узлов интерполяции берутся корни многочленов Чебышева (приведенные к отрезку [а, b], то
Чебышевской нормой функции f(x) в пространстве С[а, b] непрерывных на [а, b] функций называется
Для достаточно гладких функций f(x) при специальном расположении узлов интерполяции последовательность интерполяционных многочленов Лагранжа Ln(x) (построенных по точным значениям функции f(x) сходится к f(x) по норме Чебышева; другими словами, имеет место равномерная сходимость последовательности Ln(x) (к f x)).
Обобщением установленного факта для непрерывных (не обязательно дифференцируемых) функций и произвольных (не обязательно интерполяционных) многочленов является широко известная в математическом анализе теорема.
Теорема 18.1 (Вейерштрасса). Для любой непрерывной на [а, b] функции f(x) найдется такой многочлен Qn(x), что
Как следует из теоремы Вейерштрасса, если отказаться от того, чтобы аппроксимирующий функцию многочлен Qn(x) степени п был интерполяционным, от f(x) достаточно потребовать непрерывность на [а, b], чтобы за счет повышения степени многочдена при надлежащем подборе его коэффициентов величина чебышевской нормы могла быть сделанной сколь угодно малой, иначе, чтобы качество аппроксимации функции f(x) многочленом Qn(x) на отрезке [а, b] было сколь угодно хорошим в смысле чебышевской метрики.
Если степень п аппроксимирующего f(x) многочлена Qn(x) зафиксировать и распоряжаться только его коэффициентами, то в общем случае величину нельзя сделать сколь угодно малой. Однако доказано, что для любой функции f(x)С[а, b] существует единственный такой многочлен Qn(x), который из всех многочленов степени и наилучшим образом аппроксимирует на [а, b] функцию f(x), минимизируя максимальное расстояние между f(x) и Qn(x). Этот многочлен, т.е. многочлен , такой, что
называется многочленом наилучшего равномерного приближения для f(x) на [а, b] или ее чебышевским приближением.
Метод наименьших квадратов (также смотри шпору 42):
Из каких-либо соображений (аналитических, графических или иных) аппроксимирующая f(x) функция берется из определенного т-параметрического семейства функций, и ее параметры подбираются так, чтобы сумма квадратов отклонений вычисляемых значений (xi) от заданных приближенных значений уi была минимальной. Такая функция (т.е. при таком оптимальном наборе параметров) будет наилучшей аппроксимацией f(x) среди функций выбранного семейства в смысле метода наименьших квадратов. Ясно, что число данных приближенных значений yi в таблице должно быть не меньшим, чем число параметров в подбираемой зависимости (x); как правило, считается, что п >> т.
Итак, согласно МНК, задаем семейство у = (х,а1,а2,...,ат) и ищем значения параметров а1, а2 ,... , ат, решая экстремальную задачу
Oоптимальный набор параметров может быть найден из системы уравнений
(18.4)
представляющей необходимые условия экстремума функции Ф(а1, а2,..., ат), в силу ее специфики, являющиеся и достаточными условиями ее минимума.
Если функция (х,а1,а2,...,ат) есть линейная функция относительно своих параметров а1, а2,..., ат, то система (18.4) тоже будет линейной; в общем случае (18.4) — нелинейная система, что влечет за собой определенные трудности при ее решении. Спасительным в последней ситуации является тот факт, что обычно при задании семейств функций, аппроксимирующих реальные зависимости, число параметров берется небольшим 2 или 3, причем какие-то из этих параметров могут входить линейным образом.
В зависимости от характера табличных данных, изучаемого с помощью их изображения в соответствующей системе координат или с помощью некоторых прикидочных расчетов, при обработке результатов экспериментов часто используют те или иные из следующих двухпараметрических семейств функций:
y = ax + b, y = a + blnx (y = a + blgx),
y = axb, y = aebx (y = a10bx),
реже применяются трехпараметрические семейства
у = ах2 +bх + с, у = ахb +с, y = aebx+c ;
при изучении периодических явлений применяют тригонометрические функции.
Заметим, что вместо того, чтобы решать нелинейные системы, получающиеся из (18.4) при поиске параметров конкретных семейств функций, когда эти параметры входят туда нелинейным образом, можно попытаться сначала линеаризовать подбираемую зависимость.
26.Метод наименьших квадратов.
Предположим, что между независимой переменной х и зависимой переменной у имеется некая неизвестная функциональная связь у = f(x). Эта функция задана таблично
|
x |
x0 |
x1 |
… |
xn |
|
y |
y0 |
y1 |
… |
yn |
приближенных значений уi = f(xi), получаемых в ходе наблюдений или экспериментов. Требуется дать приближенное аналитическое описание этой связи, т.е. подобрать функцию (x) такую, которая аппроксимировала бы на отрезке [x0, xn] заданную отдельными приближенными значениями уi функцию f(xi).
Для решения этой задачи заведомо неудачным является интерполяционный подход хотя бы потому, что функция (x) такая, что (xi) = уi (при всех i{0,1,..., …, n}), будет мало похожа на искомую f(xi), поскольку в ней отразятся все ошибки экспериментальных данных. Уже это заставляет отказаться от идеи интерполяции и находить функцию (x) такую, чтобы она хорошо отражала «в среднем» зависимость между х и у.
Конкретнее, из каких-либо соображений (аналитических, графических или иных) аппроксимирующая f(x) функция берется из определенного т-параметрического семейства функций, и ее параметры подбираются так, чтобы сумма квадратов отклонений вычисляемых значений (xi) от заданных приближенных значений уi была минимальной. Такая функция (т.е. при таком оптимальном наборе параметров) будет наилучшей аппроксимацией f(x) среди функций выбранного семейства в смысле метода наименьших квадратов. Ясно, что число данных приближенных значений yi в таблице должно быть не меньшим, чем число параметров в подбираемой зависимости (x); как правило, считается, что п >> т.
Итак, согласно МНК, задаем семейство у = (х,а1,а2,...,ат) и ищем значения параметров а1, а2 ,... , ат, решая экстремальную задачу
Oоптимальный набор параметров может быть найден из системы уравнений
(18.4)представляющей необходимые условия экстремума функции Ф(а1, а2,..., ат), в силу ее специфики, являющиеся и достаточными условиями ее минимума.
Если функция (х,а1,а2,...,ат) есть линейная функция относительно своих параметров а1, а2,..., ат, то система (18.4) тоже будет линейной; в общем случае (18.4) — нелинейная система, что влечет за собой определенные трудности при ее решении. Спасительным в последней ситуации является тот факт, что обычно при задании семейств функций, аппроксимирующих реальные зависимости, число параметров берется небольшим 2 или 3, причем какие-то из этих параметров могут входить линейным образом.
В зависимости от характера табличных данных, изучаемого с помощью их изображения в соответствующей системе координат или с помощью некоторых прикидочных расчетов (см., например, [32, 36]), при обработке результатов экспериментов часто используют те или иные из следующих двухпараметрических семейств функций:
y = ax + b, y = a + blnx (y = a + blgx), y = axb, y = aebx (y = a10bx),
реже применяются трехпараметрические семейства у = ах2 +bх + с, у = ахb +с, y = aebx+c ;
при изучении периодических явлений применяют тригонометрические функции
27. СРЕДНЕКВАДРАТИЧНЫЕ ПРИБЛИЖЕНИЯ
Будем параллельно рассматривать: а) пространство СL[а, b] непрерывных на [а, b] функций со скалярным произведением (19.1) и нормой б) пространство Rn+1[а, b] сеточных функций, определенных в точках хi [а, b] (i = 0,1, .... n), со скалярным произведением (19.2) и нормой . Можно определить расстояние м\у двумя функциями . Введением расстояний м\у функциями можно исследовать близость двух функций. будут близки, если . . называют среднеквадратичным отклонением от f и f от . Для f следует найти так, чтобы . (x) называется наилучшим среднеквадратическим приближением f(x) на заданном семействе функций.
Рассмотрим, к чему сводится процесс построения наилучших среднеквадратических приближений в одном конкретном, но достаточно общем случае, когда аппроксимирующая f(x) функция (x) представляет собой линейную комбинацию нескольких других, вообще говоря, более простых (базисных) функций. Пусть — некоторая заданная на [а, b] система линейно независимых функций. Обобщенным многочленом будем называть функцию Qm (x) := с11(x) + ... + ст1(x) (19.3). Следует подобрать коэф-ты сn так, чтобы , . Т.о., надо минимизировать выражение: , которое после элементарных преобразований может быть переписано в терминах скалярных произведений (см. (19.1) и (19.2)):
(19.4)
Если функции образуют систему линейно независимых элементов пространства, то полученная симметричная линейная алгебраическая система, называемая нормальной системой МНК, имеет заведомо отличный от нуля определитель (это известный определитель матрицы Грамма, и значит, однозначно разрешима. Анализируя СЛАУ (19.4), приходим к выводу, что она упрощается в случае, когда базисные функции образуют на [а, b] ортогональную систему, т.е. (j, k)= 0 при любых k j. Следовательно, после упрощения получим:
(19.5)
В таком случае найденные коэффициенты называются коэффициентами Фурье, а обобщенный многочлен (19.3) называется обобщенным многочленом Фурье.
Еще проще построение таких многочленов, когда система — ортонормированная. Тогда ||j|| = 1, и из (19.5) следует, что cj = (j, f) при любом j {0,1,...,m}, т.е. аппроксимация функции f(x) обобщенным многочленом Фурье имеет вид
(19.6)
В частности, в математическом анализе достаточно подробно изучаются представления функций (не обязательно непрерывных) выражениями типа (19.6), в которых в качестве базисных функций используются функции ортогональной на [ , ] системы {sin kx, cos kx | k N0}. Такие представления в этом случае называются отрезками тригонометрических рядов Фурье.
28. Ортогональные многочлены и их применение для аппроксимации функций. Многочлены Лежандра Pn(x) – наиболее употребительные из классических ортогональных многочленов и единственные, для которых условие их ортогональности на отрезке [–1, 1] выполняется «в чистом виде», т.е.
Известна формула Родрига (2)
которая дает возможность в явном виде представить многочлены Лежандра.
Положив =1 и найдя по формуле (2) P1 =х, все последующие многочлены Лежандра можно получать один за другим с помощью рекуррентной формулы
При n = 1, 2, 3,... из нее имеем соответственно:
и т.д.
Сравнивая многочлены Лежандра с многочленами Чебышева, можно обнаружить, что их сходство — не только внешнее. Эти многочлены характеризует ряд одинаковых свойств. Например, как и многочлены Чебышева, многочлены Лежандра n-й степени имеют на отрезке [−1, 1] ровно n различных действительных корней. В то же время, как и многочлены Лежандра, многочлены Чебышева Тn(х) относятся к числу классических ортогональных многочленов, но ортогональность здесь понимается в более широком смысле — это ортогональность с весом. А именно, условие ортогональности многочленов Чебышева Тn(х) на отрезке [−1, 1] имеет вид
Многочлены Лагерра Ln(x) определяются требованием их ортогональности на промежутке [0, +∞) с весовой функцией , а именно условием
Многочлены Лагерра удовлетворяют рекуррентному соотношению
Так как L0 = 1, L1= −х +1. Отсюда легко получить несколько первых многочленов Лагерра: L2(x) = L3(x) = L4=
Многочлены Эрмита Нn(х) ортогональны на всей числовой оси с весом .Они удовлетворяют интегральному условию ортогональности
и рекуррентному соотношению
Так как H0 =1, H1= 2х. Отсюда можно получить несколько первых многочленов Эрмита: H2(x) = H3(x)= H4(x) =и т.д.
29. Применение интерполирования к вычислению производных
Источником формул численного дифф-ния явл. полиномиальная интерполяция. Зная в точках xi :=xi+ih (i = 0,1,..., п ) при некотором h > 0 значения yi:=f(xi) функции f(x), можно найти конечные разности kyi и записать для нее, например, первый интерполяционный многочлен Ньютона Рn(х). Дифференцируя приближенное равенство f(x)Рn(х), будем строить формулы приближенного дифф-ния разной точности в зависимости от степени n используемого интерполяционного многочлена. Через вспомогательную переменную приближенное представление функции f(x) по первой формуле Ньютона выглядит наиболее просто:
(20.1)
Отсюда получаем конечноразностную формулу численного дифференцирования
(20.2)
При использовании последнего равенства для приближенного вычисления производной функции f(x) в заданной точке из некоторой окрестности точки х0 следует найти соответствующее значение переменной и подставить его в формулу (20.2). Максимальный порядок конечных разностей в этой формуле при желании получить производную с наибольшей точностью определяется в конкретной ситуации в соответствии с приведенными в §§ 17.1, 17.2 соображениями о выборе подходящей степени интерполяционного многочлена.
Рассм. несколько частных случаев формулы (20.2), фиксируя степень n лежащего в ее основе интерполяционного многочлена (20.1) равной 1, 2, 3. Этим значениям n отвечают соответственно 1, 2, 3 первых слагаемых в формуле (20.2). Т. о., имеем:
на основе линейной интерполяции
(20.3)
на основе квадратичной интерполяции
(20.4)
на основе кубической интерполяции
(20.5)
и т.д. (при некоторых > 0, определяющих промежуток экстраполяции приемлемого качества соответствующей интерполяционной формулой).
Для дальнейшего особый интерес представляют частные случаи формул (20.3-5), связывающие приближ. значение производной f(x) в узлах х0, х1, ... с узловыми значениями самой функции. Учитывая, что точкам х0, х1, х2, х3, соответствуют значения q = 0, 1, 2, 3, и раскрывая конечные разности через значения yi (i = 0, 1, 2, 3), имеем:
при n = 1 из (20.3)
(20.6)
(20.7)
при n = 2 из (20.4)
(20.8)
(20.9)
(20.10)
при и = 3 из (20.5)
Повторное дифференцирование приближенного равенства (20.2), т.е. взятие производной по x от правой части формулы (20.2) с учетом приводит к конечноразностной формуле вычисления второй производной
(20.11)
из кот. таким же образом следует приближенная формула для третьей производной
Наиболее важной в приложениях является простейшая аппроксимация второй производной с помощью первого слагаемого, получающаяся из (20.11). В частности, в точке x1 имеем приближенное равенство
(20.12)
кот. вместе с формулами (20.6) – (20.10) широко используется при построении конечноразностных методов решения обыкновенных диффуров второго порядка.
30.Погрешность формул приближенного дифференцирования.
(20.6); (20.7); (20.8); (20.9); (20.10); (20.12)
Чтобы получить представление о точности простейших аппроксимаций значений производных в узловых точках, определяемых формулами (20.6) – (20.10), (20.12), будем предполагать, что данная функция f(x) обладает достаточной для выведения остаточных членов гладкостью. Знание структуры приближенных выражений для производных, полученных из интерполяционных соображений, позволяет без особого труда (по крайней мере, для симметричных аппроксимаций) вывести формулы их остаточных членов, манипулируя разложениями f(x) по формуле Тейлора подходящих порядков.
Простейшая несимметричная аппроксимация f(xi) (формулы первого порядка точности). Запишем представление функции f(x) по формуле Тейлора в окрестности точки xi: Выразив отсюда f(x), имеем (20.13) Первый член правой части этого равенства – разностное отношение, аппроксимирующее производную вблизи xi, а второй – остаточный член, характеризующий точность такой аппроксимации. При фиксировании в (20.13) х = xi одновременно зафиксируется и неизвестная точка ; таким образом, приходим к формуле левой аппроксимации f(xi) с остаточным членом: (20.14)
Аналогично при х = xi+1 из (20.13) получаем формулу правой аппроксимации f(xi) с остаточным членом: (20.15). В приближенных равенствах (20.16) при i=0, (20.17) при i=1 узнаём выведенные ранее формулы (20.6), (20.7), а остаточные члены в (20.14), (20.15) говорят о том, что, пользуясь аппроксимациями (20.16), (20.17), мы совершаем ошибку O(h), т.е. эти формулы имеют первый порядок точности. Определенную инфо об ошибках левой и правой аппроксимаций первого порядка дает знание знаков остаточных членов.
Простейшая симметричная аппроксимация f(xi) (формула второго порядка точности). Из разложения при х=xi+1 и х=xi1 имеем соответственно и Выполнив почленное вычитание двух последних равенств, получаем откуда с помощью теоремы о среднем, примененной к сумме третьих производных в квадратных скобках, приходим к формуле симметричной аппроксимации f(xi) с остаточным членом:(20.18) где i – некоторая точка интервала (xi1, xi+1). «Основная» часть формулы (20.18) — (20.19), a вид ее остаточного члена означает, что аппроксимация имеет второй порядок точности относительно шага h.
Простейшие аппроксимации второй производной. Из представления имеем откуда почленным сложением получаем
Выражая из последнего равенства , приходим к формуле симметричной аппроксимации с остаточным членом: (20.20)
Остаточный член этой формулы некоторым характеризует приближенное равенство (20.21) как аппроксимацию второй производной в точке хi, второго порядка точности, т.е. с погрешностью O(h2).
То же отношение (20.21), используемое в качестве несимметричной аппроксимации второй производной функции f(x), т.е. для вычисления приближенных значений и дает лишь первый порядок точности.
31. Квадратурные формулы и связанные с ними задачи. Интерполяционные квадратурные формулы прямоугольников, трапеций.
Пусть требуется найти значение I интеграла Римана для некоторой заданной на отрезке [а, b] функции f(х).
Хорошо известно, что для функций, допускающих на промежутке [а, b] конечное число точек разрыва первого рода, такое значение существует, единственно и может быть формально получено по определению: (21.1) где — произвольная упорядоче нная по возрастанию система точек отрезка [а, b] такая, что а i — произвольная точка элементарного промежутка [xi1, xi].
Простые квадратурные формулы можно вывести непосредственно из определения интеграла, т.е. из представления (21.1). Зафиксировав там некоторое n > 1, будем иметь приближенное равенство, которое назовем общей формулой прямоугольников (площадь криволинейной трапеции приближенно заменяется площадью ступенчатой фигуры, составленной из прямоугольников, основаниями которых служат отрезки [xi1, xi], а высотами — ординаты f(i).
Чтобы получить конструктивное правило приближенного вычисления интеграла, воспользуемся свободой расположения точек xi разбивающих промежуток интегрирования [а, b] на элементарные отрезки [xi1, xi], и свободой выбора точек i,- на этих отрезках.
Условимся в дальнейшем в данном прараграфе пользоваться только равномерным разбиением отрезка [а, b] на n частей точками xi с шагом , полагая х0=а, хi = хi1 +h (i = 1,2, ...,n1), хn=b. (21.2)
При таком разбиении формула (21.1) приобретает вид
Теперь дело за фиксированием точек i, на элементарных отрезках [xi1, xi]. Рассмотрим три случая.
1) Положим i = xi1. Тогда получаем (21.4)
2) Пусть в (21.3) i = xi. В таком случае имеем (21.5)
Формулы (21.4) и (21.5) называются квадратурными формулами левых и правых прямоугольников соответственно. Совершенно очевидно, что и в совокупности дают двусторонние приближения к значению I интеграла от монотонной функции.
Можно рассчитывать на большую точность получения значения интеграла, если взять точку i, посередине между точками xi1 и xi. Отсюда приходим к третьему случаю.
3) Фиксируем В результате имеем квадратурную формулу средних прямоугольников или, иначе, формулу средней точки (21.6)
Пусть функция f(x) дважды непрерывно дифференцируема. Рассмотрим сначала вычисление интеграла (21.7) при некотором достаточно малом h > 0. По формуле Тейлора можно записать , где х — произвольная, a — некоторая фиксированная точки интервала (h/2, h/2). Подставляя полученное выражение для f(x) в (21.7), имеем:
(21.8) где 0 (h/2, h/2) — некоторая точка, вообще говоря, несовпадающая с точкой (поскольку значение изменяется с измененнием h)
Слагаемое hf(0) в выражении (21.8) интеграла I0, очевидно, можно расценивать как приближение к I0 по формуле прямоугольников в случае всего одного элементарного отрезка. Второе же слагаемое, т.е. есть остаточный член (локальная погрешность) простейшей формулы прямоугольников I0=hf(0).
В результате подстановки этого выражения в каждое слагаемое правой части равенство приходим к окончательному виду остаточного члена (глобальной погрешности) квадратурной формулы прямоугольников:
(21.9)
Как видно из формулы (21.9), при увеличении числа n элементарных отрезков, на которые разбивается промежуток интегрирования [а, b], ошибка численного интегрирования по формуле средней точки (21.6) убывает пропорционально квадрату шага h. Нетрудно убедиться, что погрешность численного интегрирования непрерывно дифференцируемой функции по формулам левых и правых прямоугольников (21.4), (21.5) убывает лишь по линейному закону.
32.Простейшие квадратурные формулы Ньютона-Котeса. Квадратурная формула Симпсона. Оценки точности квадратурных формул.
Подстановка в интеграл вместо функции f(x) ее интерполяционного многочлена Лагранжа той или иной степени n приводит к семейству квадратурных формул, называемых формулами Ньютона-Котеса.
функция может быть единственным образом представлена в виде
f(x) = Ln(x) + Rn(x), где Ln(x) — интерполяционный многочлен Лагранжа, Rn(x) — остаточный член.
Если ситема узловинтерполирования совпадает с точками разбиения (21.2) отрезка [а, b] с шагом h, то замена переменной x = x0+qh трансформирует многочлен Лагранжа к виду
Для того чтобы использовать такое выражение Ln(x) вместо f(x) в , нужно изменить границы интегрирования (значению х = а соответствует значение q = 0, а х = b — значение q = n и учесть, что dх = hdq. Таким образом, получаем
Это равенство, переписанное в виде и есть квадратурная формула Ньютона—Котеса, где — коэффициенты Котеса.
1) Пусть n = 1, т.е. имеется всего две точки x0 и x1 = x0 +h, в которых известны значения функции (y0 = f(x0) и yi = f(x1). Этим точкам соответствуют значения 0 и 1 переменной q. Следовательно,
Получена простейшая квадратурная формула трапеций, к которой легко прийти и из геометрических соображений.
Остаточный член этой формулы найдем интегрированием остаточного члена R1(x) формулы линейной интерполяции, преобразованного к виду
Учитывая, что dх=hdq имеем:
В найденном выражении остаточного члена 1 = (x0, x1)— некоторая точка, отвечающая интегральной теореме о среднем, которая здесь также применима, в силу неположительности функции q2 - q на отрезке [0,1].
2) Положим в (21.12) n = 2, т.е. проинтерполируем функцию /(*) по трем точкам: x0, x1 = x0 +h и x2 = x0 + 2h. Тогда
Полученное приближенное равенство назовем простейшей формулой Симпсона.
Аналогично как и в п. 1) можно найти остаточный член простейшей формулы Симпсона
Фиксируя степень интерполяционного многочлена (21.12) равной 3, 4, 5, …, получим частные формулы Ньютона-Котеса, подобно полученным выше простейшим формулам трапеций и Симпсона.
33. Принцип Рунге практического оценивания погрешностей.
Пусть к приближенному вычислению значения f данного интеграла применяется некая квадратурная формула р-го порядка точности IР из семейства составных формул Ньютона-Котеса. При условии непрерывности p-й производной подынтегральной функции это означает существование такой постоянной С, что
I = Ip(h) + Chp. (21.13)
При уменьшении вдвое шага h численного интегрирования по той же формуле р-го порядка можно записать такое же равенство, но с другой постоянной C1 :
I =Ip(h/2) + C1(h/2)p. (21.14)
Считая, что при малом h постоянные С и С1 близки, из (21.13) и (21.14) имеем
Ip(h) + Chp Ip(h/2) + C(h/2)p.
и, следовательно,
Подставив это значение С1 в (21.14), приходим к выражению
(21.15)
К приближенному равенству (21.15) можно относиться двояко. С одной стороны, переписав его в виде
(21.16)
получаем возможность хотя бы грубо контролировать точность численного интегрирования на основе двойного счета (с шагом h и с шагом h/2). В этом и заключается широко применяемый принцип Рунге практического оценивания погрешностей. Его применение считается правомочным, если выполняется неравенство
[[[ этого нет в конспекте: С другой стороны, второе слагаемое в формуле (21.15) позволяет уточнить «дешевым» способом приближенное значение IР (h/2) интеграла I. Заметим, что при р = 2 формула (21.13) есть формула трапеций
I = IT(h) + Ch2,
и выражение (21.16) в этом случае является поправкой Ричардсона RT (h/2); а правая часть формулы (21.15) есть значение Ic(h/2), соответствующее формуле Симпсона. Естественно назвать обобщенной поправкой Ричардсона ( или как в конспекте просто ОШИБКА) получаемую с помощью двойного счета величину ]]]
Ошибка:
34. Методы решения задачи Коши. Метод Эйлера — геометрический подход к построению.
Учитывая ключевую позицию, которую занимает метод Эйлера в теории численных методов ОДУ, рассмотрим несколько способов его вывода. При этом будем считать, что вычисления проводятся с расчетным шагом расчетными точками (узлами) служат точки хi = х0 + ih
(i = 0,1,.,., n) промежутка [х0, b] и целью является построение таблицы
|
x |
x0 |
x1 |
… |
xn=b |
|
|
y |
y0 |
y1 |
… |
yn=y(b) |
приближенных значений yi решения у = у(х) задачи (1)
y(x0) = y0 (2) в расчетных точках xi.
Геометрический способ. Пользуясь тем, что в точке x0 известно и значение решения y(x0) = y0, и значение его производной , можно записать уравнение касательной к графику искомой функции у = у(х) в точке (х0;у0): (3)При достаточно малом шаге h ордината (4)
этой касательной, полученная подстановкой в правую часть (3) значения х1 = х0 + h по теореме о непрерывной зависимости решений от начальных данных должна мало отличаться от y(x1) ординаты решения y(x) задачи (1)-(2). Следовательно, точка (x1,y1) может быть приближенно принята за новую начальную точку. Через эту точку снова проведем прямую , которая уже приближенно отражает поведение касательной к решению у = у(х) в точке. Подставляя сюда + h), получим приближение значения у(x2) значением и т.д. В итоге этого процесса, определяемого формулой (5)
и называемого методом Эйлера, график решения у = у(х) данной задачи Коши (1) – (2) приближенно представляется ломаной, составленной из отрезков приближенных касательных, откуда происходит другое название метода (5) —метод ломаных.
35. Применение формулы Тейлора к построению метода Эйлера.
Применение формулы Тейлора. Линеаризуя решение в окрестности начальной точки по формуле Тейлора, имеем Отсюда при x = x1 получаем Сравнивая это равенство с формулой (4), видим, что ее остаточный член имеет вид (6) где — некоторая точка интервала (x0, х1). Остаточный член (6) характеризует локальную (или, иначе, шаговую) ошибку метода Эйлера, т.е. ошибку, совершаемую на одном шаге. Очевидно, что от шага к шагу, т.е. при многократном применении формулы (5), возможно накопление ошибок. За n шагов, т.е. в точке b, образуется глобальная ошибка. Порядок глобальной ошибки (относительно шага h) на единицу ниже, чем порядок локальной ошибки, а порядком глобальной ошибки и определяется порядок соответствующего численного процесса решения задачи Коши. Таким образом, локальная ошибка метода Эйлера, согласно (6), есть O(h2), глобальная — O(h), т.е. метод Эйлера относится к методам первого порядка.
36. Квадратурный способ построения метода Эйлера
Квадратурный способ. Начальную задачу для ОДУ (1)-(2) можно заменить эквивалентным интегральным уравнением (7) При x = x1 из него получится равенство (8) Применение к интегралу в правой части равенства простейшей (одноточечной) формулы левых прямоугольников дает приближенную формулу
правая часть которой, очевидно, совпадает с выражением (4) для подсчета значения y1. Для получения формулы (5) следует применить метод левых прямоугольников к интегралу (9)которая получается из формулы (8) заменой точки (x0, y0) точкой (xi, yi), а точки (x1, y1) точкой (xi+1, yi+1).
37.Одношаговые методы типа Рунге-Кутта
Недостатком методов более высоких порядков, основанных на пошаговом представлении решения у(х) по формуле Тейлора и последовательном дифференцировании уравнения для получения тейлоровых коэффициентов, является необходимость вычисления на каждом шаге частных производных функции f(х, у).
Идея построения явных методов Рунге-Кутты р-ro порядка заключается в получении приближений к значениям по формуле вида (22.12)
где j(х, у, h) — некоторая функция, приближающая отрезок ряда Тейлора до р-гo порядка и не содержащая частных производных функции f(x, у).
Так, полагая в (22.12) j(х, у, h) = f(x, у), приходим к методу Эйлера (22.5), т.е. метод Эйлера можно считать простейшим примером методов Рунге-Кутты, соответствующим случаю р = 1.
Для построения методов Рунге—Кутты порядка, выше первого, функцию j(х, у, h) берут многопараметрической и подбирают ее параметры сравнением выражения (22.12) с многочленом Тейлора для у(х) соответствующей желаемому порядку степени.
Рассмотрим случай р = 2. Возьмем функцию j в (22.12) следующей структуры:
Ее параметры с1, с2, а и b будем подбирать так, чтобы записанная, согласно (22.12), формула
(22.13)
определяла метод второго порядка, т.е. чтобы максимальная локальная ошибка составляла величину O(h3).
Разложим функцию двух переменных f(x + ah, y + bhf(x, у)) по формуле Тейлора, ограничиваясь линейными членами:
Ее подстановка в (22.13) дает
Сравнение последнего выражения с тейлоровским квадратичным представлением решения у(х):
с точностью до O(h3) от параметров нужно потребовать выполнение следующей совокупности условий:
(22.14)
Полученная система условий содержит три уравнения относительно четырех параметров метода. Это говорит о наличии одного свободного параметра. Положим с2 =a (¹ 0). Тогда из (22.14) имеем:
В результате подстановки этих значений параметров в формулу (22.13) приходим к однопараметрическому семейству методов Рунге—Кутты второго порядка.
(22.15)
Выделим из семейства методов (7.30) два наиболее простых и естественных частных случая:
при а = 1/2 получаем формулу , которая называют методом Хойна;
при a = 1 из (7.30) выводим новый простой метод ,который назовем методом средней точки.
38.Методы Рунге-Кутты произвольного и четвертого порядков
Любой метод из семейства методов Рунге-Кутты второго порядка (7.30) реализуют по следующей схеме. На каждом шаге, т.е. при каждом = 0, 1, 2,..., вычисляют значения функции: а затем находят шаговую поправку
прибавление которой к результату предыдущего шага дает приближенное значение решения у(х) в точке :
Метод такой структуры называют двухэтапным по количеству вычислений значений функции — правой части уравнения— на одном шаге.
Анализ устройства методов Рунге-Кутты второго порядка позволяет представить, в какой форме следует конструировать явный метод Рунге-Кутты произвольного порядка. По аналогии с предыдущим для семейства методов Рунге-Кутты р-го порядка используется запись, состоящая из следующей совокупности формул:
где k = 2,3,..., p (для р-этапного метода). Многочисленные параметры фигурирующие в формулах (7.32), подбираются так, чтобы получаемое методом (7.32) значение совпадало со значением разложения у(хi+1) по формуле Тейлора с погрешностью (без учета погрешностей, совершаемых на предыд. шагах).
Наиболее употребительным частным случаем семейства методов является метод Рунге-Кутты четвертого порядка, относящийся к четырехэтапным и определяемый следующими расчетными формулами*:
Δ
39. Линейные многошаговые методы
Пусть требуется найти решение u(t) на отрезке [t0> Т] задачи Коши Предположим, что построены приближенные значения решения u(t) и его первой производной u'(t) в моменты времени ,т. е.
Общий вид разностной схемы рассматриваемых здесь многошаговых методов имеет вид где коэффициенты (их всего 2p + 3), которые должны быть определены при получении конкретного многошагового метода, h - шаг интегрирования.
Значения этих 2p + 3 коэффициентов выбирают так, что если решение u(t) является полиномом степени n, то разностная схема многошагового метода дает точное значение, т. е. yk+1 =u(tk+1). Поскольку полином степени n u(t) = a0 +a1t + a2t2 +--- + antn имеет n + 1 параметр, то разностная схема должна иметь по крайней мере п + 1 коэффициент. В большинстве практических многошаговых методов 2р + 3>п + 1 и лишние коэффициенты могут быть выбраны произвольно.
Получим соотношения, которым должны удовлетворять все 2р + 3 коэффициента разностной схемы в предположении, что метод дает точное решение для задачи Коши, точным решением которой является полином степени п . 1. п = 0 , u(t) = а0 . Класс задач с таким решением задается уравнением Поэтому Подставив эти значения в разностную схему, получим первое условие, которому должны удовлетворять коэффициенты : 2. п = 1, u(t) = a0 +a1t. Класс задач с таким решением задается в виде Для удобства выберем tk =0. Тогда tk+1 = h , tk1 = h, tk2 = 2h, tki = hi. В этом случае Подставим их в разностную схему: Если учесть первое условие, то это соотношение преобразуется к виду Наконец, разделив левую и правую части на a1h, получим условие корректности для полиномиальных решений первой степени: 3. n = 2, u(t) = a0 +a1t+ a2t2. Класс задач с таким решением Полагаем, как и прежде, tk = 0, tki = hi и находим, что Перепишем с учетом этих соотношений разностную схему многошагового метода Разделим левую и правую части этого соотношения на . Условие корректности для полиномиальных решений второй степени примет следующий вид: 4. Общий случай: u(t) = a0 +a1t+ a2t2+ …+antn Класс задач с таким решением решения до степени п включительно имеют одинаковую форму, а именно: Этим соотношениям должны удовлетворять все 2p + 3 коэффициента разностной схемы линейного многошагового метода.
40. Явные многошаговые методы Адамса.
Этот класс методов легко получить из общей разностной схемы линейных многошаговых методов и условий корректного выбора коэффициентов Полагая р = n 1, а1= а2 =… = ап = 0, b1 = 0, получим разностную схему явного многошагового метода порядка n где коэф-ты,(их (n+1)),однозначно опред-ся из (n+1)условия корректности. Из первого условия корректности находим а0 =1.
Обозначим bi = . Эта означает, что значение коэфф-та bj зависит от порядка метода. Чтобы определить подставим значения ,= 1, аi=0, в оставшиеся п условий корректности: Получили систему из n лин-х алгебр-х уравнений относительно
.
Запишем явные методы Адамса различных порядков.
Разност. схема явн. метода Адамса 1го порядка совпадает с разност. схемой явн. метода Эйлера.
2. n = 2. Соотношения для являются системой из 2х лин-х алгебр-х уравнений Следовательно, (2) = 1/2,(2) = 3/2, и разност. схема явн. метода Адамса 2го порядка выглядит
3. = 3. Относительно (3), = 0,2, запишем систему лин-х алгебр-х уравнений 3го порядка из которой следует, что (3) = 5/12, (3) =16/12. (3) = 23/12 . В результате разностная схема явного метода Адамса третьего порядка записывается след. образом: Этот процесс записи разностных схем явных многошаговых методов более высоких порядков точности можно продолжить и далее. Не останавливаясь на выводе, приведем выражение для локальной погрешности явных многошаговых методов Адамса порядка n : где – постоян. коэф-т, величина которого зависит от метода, в частности =1/2, С2 =5/12, С3 =3/8;
-значение (+1)-й производной функции u(t) в точке
41.Неявные многошаговые методы Адамса.
Неявный метод Адамса n -го порядка получается из общей разностной схемы линейных многошаговых методов при условии: р = n 2, а1= а2 =… = аn2 = 0,
т. е.
Здесь a0,b1, b0 b1, bn2 их всего n+1, определяются из условий корректности полиномиальных решений n-го порядка.
Как и в явных методах Адамса, а0 =1 . Обозначим bi = i(n), подчеркивая тем самым зависимость значений коэффициентов bi от порядка метода. Подставим выбранные значения параметров в условие корректности:
Относительно неизвестных коэффициентов i (n), получим систему из п линейных алгебраических уравнений:
Решение этой системы однозначно определяет коэффициенты неявного метода Адамса n -го порядка.
Запишем неявные методы Адамса первого, второго и третьего порядков.
1. n = 1. Коэффициент 1 (1) = 1 и
Разностная схема неявного метода Адамса первого порядка совпадает с неявным методом Эйлера.
2. n = 2 . В этом случае: и искомые коэффициенты имеют следующие значения: .
Неявный метод Адамса второго порядка: как легко видеть, является методом трапеций.
3. n = 3. Коэффициенты вычисляются в этом случае из системы: Ее решение .
Неявный метод Адамса 3-го порядка принимает вид: И т.д.
Приведем теперь без вывода соотношение для расчета локальной погрешности неявного метода Адамса n-го порядка: , где значения коэффициента Сn для рассмотренных выше разностных схем соответственно равны: С1 = 1/2, С2 = 1/12, С3 = 1/24.
Сравнивая явные и неявные методы Адамса одного порядка между собой, можно отметить, что методы имеют одинаковую по порядку, но разную по знаку локальную погрешность.
42. Условие устойчивости линейного многошагового метода.
Будем полагать теперь, что для решения тестового уравнения используется многошаговый метод. В этом случае
Полагая σ=hλ, преобразуем это выражение следующим образом:
(1+ σb-1)yk+1-(a0- σb0)yk-(a1- σb1)yk-1-…-(ap- σbp)yk-p=0.
Сделаем подстановку: yk=zk.
Тогда разностное уравнение приводится к виду
(1+ σb-1)zk+1-(a0- σb0) zk-(a1- σb1) zk-1-…-(ap- σbp) zk-p=0.
Или: zk-p ((1+ σb-1)zp+1-(a0- σb0) zp-(a1- σb1) zp-1-…-(ap- σbp)).
Каждый корень полинома
P(z)=(1+ σb-1)zp+1-(a0- σb0) zp-(a1- σb1) zp-1-…-(ap- σbp)
Их всего p+1 для каждого фиксированного σ, порождает частное решение yk=z разностного уравнения (25.1) многошагового метода. Так как разностное уравнение многошагового метода линейно и для него справедлив принцип суперпозиции. Предполагая, что все корни zk, p=1,…,p+1 различны, общее решение можно представить в виде
Yk=c1z1k+ c2z2k+…+ cp+1zp+1k ,
Где с – постоянные коэффициенты. Если полиномиальное уравнение P(z)=0 содержит кратный корень zi, кратности mi, то соответствующий член в общем решении будет таким
(ci0+ ci1k+ …+ ci(mi-1)kmi-1)zk.
Отсюда следует
Таким образом, многошаговый метод является численно устойчивым для тех значений σ=hλ, для которых (p+1) корень полиномиального уравнения P(z)=0 лежит внутри единичной окружности (модуль)z=1.
Множество всех значений σ=hλ, для которых многошаговый метод является численно устойчивым, называют областью устойчивости метода.