part2 / текст(all)

Обозначим среднюю точку отрезка a, b

через

i

xi 1 xi

.

(4.6)

2

Запишем для функции

f (x)

на каждом из отрезков xi 1 , xi

фор-

мулу Тейлора с остаточным членом в форме Лагранжа

f (x) f (i )

(x i ) f (i )

(x )2

f (i );

i

2!

(4.7)

i xi 1 , xi

Подставим в правую часть соотношения (4.5) вместо

f (x) ее пред-

ставление (4.7)

b

N

x1

x i f i

x i 2

f (x)dx

f (i )

2!

f i dx

a

i 1 xi 1

.

(4.8)

N

xi

xi

xi

x i 2

i

f (i ) dx f i

x i dx

2

f

dx

xi 1

xi 1

xi 1

i 1

Используя для вычисления

xi

x i 2

f i dx

вторую теорему о

2

xi 1

xi

i

x

0 , получим,

среднем значении функции и, учитывая,

что

dx

xi 1

что

b

N

h3 N

i

i

f x dx h

f

f

;

24 i 1

a

i 1

(4.9)

i xi 1, xi .

В силу непрерывности				f x	существует такая точка a, b что
N			i Nf
f				.				(4.10)
i 1
Используя (4.10), получаем
		b		N		h	3
		f (x)dx h f i						Nf
		a		i 1	24

31

	или, так как h		b a	,
			N
	b			N		1				b a h2 f .
	f x dx (b a)								f i
	a			i 1		N				24
	Приближенное равенство
	b			N			1
	f	x dx (b a)							f i J Nnp
							N
	a			i 1

называется квадратурной формулой прямоугольников,

узлами		a, b и коэффициентами A						1	.
	i
					i			N
Величина
		b		b a
		np				2
		RN ( f ) f x dx J N			h
				24			f ( )
		a

(4.11)

(4.12)

определяемой

(4.13)

является остаточным членом формулы прямоугольников (4.12).

Оценка остаточной погрешности формулы прямоугольников может быть записана в виде

			R ( f )					b a			h2 M	2		,	(4.14)
			N				24						1
где	M 2	max			f (x)					.
		a,b

Выражения для остаточного члена (4.13) и остаточной погрешности (4.14) показывают, что формула прямоугольников (4.12) является точной для любой линейной функции, так как вторая производная такой функции равна нулю, и, следовательно, 1 0 .

Оценим вычислительную погрешность 2 формулы прямоугольников, которая возникает за счет приближенного вычисления значений

функции f (x) в узлах i .
Пусть, например, значения f (i ) в формуле (4.12)					вычислены с
одинаковой абсолютной погрешностью * , тогда
				N
2	J N	J	N	(b a) * b a * .	(4.15)
				i 1

32

4.2. Метод трапеций

Предположим, что

f (x) C2 a;b .

Разделим отрезок a; b

на N

равных частей, тогда

b

N

xi

f (x)dx

f (x)dx .

(4.16)

a

i 1 xi 1

b; h

b a

где x a ih;i 0, N 1; x

N

i

N

Заменим функцию

f (x)

на каждом из отрезков xi 1 , xi

первой

интерполяционной формулой Ньютона первой степени

f (x) f (x

)

x xi 1

( f (x ) f (x

))

i 1

h

i

i 1

,

(4.17)

f (i )

(x x

)(x x )

2!

i 1

i

Подставляя формулу (4.17) в правую часть (4.16), интегрируя и используя вторую теорему о среднем значении функции, получим

b

N

f (x

) f (x )

h3 N

f (x)dx h

i 1

i

f

(4.18)

(i );i xi 1; xi ,

a

i 1

2

12 i 1

В силу (4.10) получаем:

b

f (x ) f (x

)

N

h2

N

0

f (xi )

f (x)dx h

(4.19)

2

12

b a f () .

a

i 1

Приближенное равенство

b	b a		f (x ) f (x		)	N
				N
J f (x)dx			0			f (xi )		JNTP
	N		2
a						i 1

называется формулой трапеции. Величина

RN ( f ) J JNTP 12h2 b a f ( ) .

(4.20)

(4.21)

является остаточным членом формулы трапеций. Оценка остаточной погрешности формулы трапеций может быть записана в виде

R ( f )			b a	h2 M	2		(4.22)
N		12				1

33

Формула трапеций, как и формула прямоугольников, является точной для любой линейной функции. Вычислительная погрешность формулы трапеций также равна

2 (b a) * .

(4.23)

Так как остаточные члены формул прямоугольников и трапеций (4.13) и (4.21) имеют противоположные знаки, то формулы (4.12) и (4.20) дают двухстороннее приближение для интеграла (4.1), то есть

	J Nпр J J NТР , если									f (x) 0 ,
	J NТР J J Nпр , если									f (x) 0 ,
В таком случае можно принять, что
		J пр		J ТР
	J		N		N		J ,			(4.24)
				2
тогда
					J пр	J ТР
	J J				N		N		.	(4.25)
						2

т.е. погрешность выражается через приближенные значения интегралов.

4.3. Метод Симпсона

Предположим, что		f (x) C4 a;b .						Разделим отрезок			a; b на
N 2k равных частей, тогда
	b		k 1 x2i 2
	f (x)dx					f (x)dx ,					(4.26)
	a		i 1	x2i
			b; h		b a				b a
где x a ih;i 0, 2k 1; x
i		2k				N			2k
Заменим функцию		f (x)	на каждом из отрезков x2i , x2i 2 длиной

2h по формуле Стирлинга второго порядка. Проводя рассуждения, аналогичные сделанным при выводе формуле трапеций, получим квадра-

турную формулу Симпсона

b

J f (x)dx

			a			, (4.27)
		h
				k	k 1
			f (x0 ) f (x2k ) 4 f (x2i 1 ) 2 f (x2i ) J2Ck
		3
				i 1	i 1

34

с остаточным членом

		R ( f ) J J C							h4		(b a) f IV ( ), (a;b) .				(4.28)
		N			2k			180
Оценка остаточной погрешности формулы Симпсона примет вид
				R ( f )				b a				h4 M	4	,	(4.29)
				N				180						1
где M	4	max	f IV (x)		.
		a;b
Вычислительная погрешность формулы Симпсона равна
				2	* (b a) .										(4.30)

Из выражения для остаточного члена формулы Симпсона следует, что она точна для многочленов третьей степени.

Вопросы для самопроверки

1.В каком случае используется численное интегрирование?

2.Постановка задачи численного интегрирования.

3.Какие существуют методы интегрирования функций?

4.Графическая интерпретация метода трапеций.

5.Как оценить погрешность метода трапеций?

6.Графическая интерпретация метода Симпсона.

7.Как оценить погрешность метода Симпсона?

8.Графическая интерпретация метода прямоугольников.

9.Как оценить погрешность метода прямоугольников?

10.Чем отличаются формулы метода трапеций и метода Симпсона?

11.Чем отличается вычисление погрешности метода трапеций и

Симпсона?

12. Каковы преимущества формулы парабол по сравнению с формулой трапеций и следствием чего являются эти преимущества?

13. Как влияет на точность численного интегрирования величина шага?

14. Можно ли добиться неограниченного уменьшения погрешности интегрирования путем последовательного уменьшения шага?

35

Тема 5. РЕШЕНИЕ ОБЫКНОВЕННЫХ ДИФФЕРЕНЦИАЛЬНЫХ УРАВНЕНИЙ (ЗАДАЧА КОШИ)

Определение. Обыкновенными дифференциальными уравнениями называются такие уравнения, которые содержат одну или несколько производных от искомой функций:

F (x, y, y ',..., y(n) ) 0 .

(5.1)

Наивысший порядок n входящей в уравнение (5.1) производной называется порядком дифференциального уравнения.

Решением дифференциального уравнения (5.1) называется всякая n раз дифференцируемая функция y (x) , которая после ее подстановки

в уравнение превращает его в тождество.

Для решения обыкновенных дифференциальных уравнений применяется метод конечных разностей. Его сущность состоит в следующем:

1.Область непрерывного изменения аргумента (например, отрезок) заменяется дискретным множеством точек – узлами. Эти узлы составляют разностную сетку.

2.Искомая функция непрерывного аргумента приближенно заменяется функцией дискретного аргумента на заданной сетке (сеточной функцией).

3.Исходное дифференциальное уравнение заменяется разностным уравнением относительно сеточной функции.

Такая замена дифференциального уравнения разностным называет-

ся его аппроксимацией на сетке (или разностной аппроксимацией).

Таким образом, решение дифференциального уравнения сводится к отысканию значений сеточной функции в узлах сетки.

Рассмотрим задачу Коши:

y f (x, y),	(5.2)
y(x0 ) y0 .	(5.3)

для определенности будем считать, что решение нужно получить для значений x x0 .

5.1. Метод Эйлера

Простейшим численным методом решения задачи Коши для обыкновенного дифференциального уравнения является метод Эйлера.

36

Заменим в левой части уравнения (5.2) производную y правой разностью. При этом значения функции y в узлах xi заменим значениями сеточной функций yi :

			yi 1 yi		f (x , y ) .						(5.4)
			hi		i	i
	Будем	считать для			простоты			узлы равноотстоящими, т.е.
hi h xi 1		xi const .			Тогда		из	равенства (5.4)			получаем
yi 1 yi h f (xi , yi ) .
	С помощью метода Эйлера значение сеточной функции										yi 1 в лю-
бом узле xi 1 вычисляется по ее значению									yi в предыдущем узле xi . В
связи с этим метод Эйлера относится к одношаговым методам.
	Для оценки погрешности на практике пользуются двойным просче-
том: с шагом h и шагом h / 2 .
	Погрешность более точного значения								yi (при шаге h / 2 ) оцени-
вают приближенно так:
					y(x ) y			y y
					i	i		i	i
где	yi – приближенное значение полученное при вычислениях с шагом
y
h , –i приближенное значение полученное с шагом h / 2 .

5.2. Метод Рунге-Кутта

Существуют и другие явные одношаговые методы. Широко распространен метод Рунге-Кутта четвертого порядка.

Если известно значение функции yi 1 в точке xi 1 , то вычисление приближенного значения yi в следующей точке xi xi 1 h производится по формулам (5.5).

K (k ) f (x

, y

),

1

i 1

i 1

K (k ) f (x

h

K (k )

, y

1

),

2

i 1

2

i 1

2

K (k ) f (x

h

, y

K2(k )

),

.

(5.5)

3

i 1

2

i 1

2

K (k ) f (x

h, y

K (k ) ),

4

i 1

i 1

3

yi 1

yi

h

K1(k ) 2K2(k ) 2K3(k ) K4(k )

6

37

Оценку погрешности метода можно получить с помощью двойного просчета по формуле

y(x ) y				1			y y
k	k		15				k k

Данный метод Рунге-Кутта требует на каждом шаге четырехкратного вычисления правой части f (x, y) уравнения (5.2).

Метод Рунге-Кутта требует большего объема вычислений по сравнению с методом Эйлера, однако это окупается повышенной точностью, что дает возможность проводить счет с большим шагом. Другими словами, для получения результатов с одинаковой точностью в методе Эйлера потребуется значительно меньший шаг, чем в методе Рунге-Кутта.

Схема Рунге-Кутта имеет ряд важных достоинств: она имеет хорошую точность; является Явной схемой (т.е. значение yi 1 вычисляется по ранее

найденным значениям за определенное число действий по определенным формулам; схема допускает расчет с переменным шагом; для начала расчета достаточно задать значение y0 . Все эти свойства очень ценны при расчетах.

5.3. Методы прогноза и коррекции (метод Адамса, метод Милна)

Суть методов: на каждом шаге вводятся два этапа, использующих многошаговые методы: 1) с помощью явного метода по известным значениям функции в предыдущих узлах находится начальное приближе-

ние y(xi ) yi в новом узле; 2) используя неявный метод, в результате

итераций находятся приближения yi 1 , yi 2 ,... .

Широко распространенным семейством многошаговых методов являются метод Адамса и метод Милна. Простейший из них совпадает с рассмотренным ранее методом Эйлера первого порядка точности. В практических расчетах чаще всего используется вариант метода Адамса, имеющий четвертый порядок точности и использующий на каждом шаге результаты предыдущих четырех. Именно его и называют обычно методом Адамса. Рассмотрим этот метод.

Пусть найдены значения yi в четырех последовательных узлах x0 , x1 , x2 , x3 . При этом имеются также вычисленные ранее значения правой части y0 , y1 , y2 , y3 . В качестве интерполяционного многочлена

можно взять многочлен Ньютона. Тогда разностная схема четвертого порядка для k-го шага метода Адамса имеет вид:

yпред y			h	(55 y'	59 y'	37 y'	9 y'	) ,	(5.6)
	k 1
k			24	k 1	k 2	k 3	k 4

38

( y' )пред

f (x , yпред )

,

(5.7)

k

k

k

yкор y

h

(9( y'

)пред 19 y'

5 y'

y'

) .

(5.8)

k 1

k

24

k 1

k 1

k 2

k 3

Точность вычислений оценивается по формуле:

yкор

y(x

)

1

yкор yпред

.

k

k

4

k

k

Разностные соотношения для k-ого шага метода Милна имеют вид:

yпред y

4h

(2y'

y'

2y'

) ,

(5.9)

k

k 4

3

k 3

k 2

k 1

( y'

)пред

f (x

, y

пред ) ,

(5.10)

k

k

k

yкор y

h

( y'

4y'

( y'

)пред ) .

(5.11)

k

k 2

3

k 2

k 1

k

Точность вычислений оценивается по формуле:

yкор y(x )				1		yкор yпред			.
k	k		29			k	k

Сравнивая методы прогноза и коррекции с методом Рунге – Кутта той же точности, отмечаем их экономичность, поскольку они требуют вычисления лишь одного значения правой части на каждом шаге. Но чтобы начать расчет методом Адамса или Милна, недостаточно знать y0 . Для начала расчета надо знать величину решения в четырех точках

x0 , x1 , x2 , x3 . Поэтому необходимо вычислить недостающие значения yk каким-либо другим методом, например, методом Рунге-Кутта. Кро-

ме того, методы прогноза и коррекции не позволяют изменить шаг в процессе расчета; этого недостатка лишены одношаговые методы.

5.4. Метод Хойна

Метод Хойна (метод Эйлера-Коши) это модифицированный метод Эйлера. Метод Хойна представляет собой простейшую предикторкорректорную схему: сначала вычисляется грубое приближение к решению по методу Эйлера, затем оно уточняется при помощи неявного метода трапеций, имеющего более высокий порядок точности.

Разностные соотношения для k-ого шага имеют вид:

yk	yk 1	h	f (xk 1 , yk 1 ) f (xk , yk 1	h f (xk 1	, yk 1 ))	.
			2

39

Оценка погрешности в точке xk , полученная с помощью двойного

					1			y y
	пересчета, имеет вид:	y(x ) y								.
		k	k	3				k	k

Как правило, при достаточно малом h итерации быстро сходятся. Если после трех-четырех итераций не произошло совпадение нужного числа десятичных знаков, то следует уменьшить шаг расчета h.

Вопросы для самопроверки

1.Что значит – решить задачу Коши для дифференциальных уравнений первого порядка?

2.Графическая интерпретация численного решения дифференциального уравнения.

3.Какие существуют методы решения дифференциального уравнения в зависимости от формы представления решения?

4.В чем заключается суть метода Эйлера?

5.В чем заключается суть метода Рунге-Кутты?

6.В чем заключается суть методов прогноза и коррекции?

7.Как вычислить погрешность по заданной формуле, используя метод двойного пересчета?

40