Учебное пособие_Численные методы -26 -дек_2014 (97-2003)

Останавливаясь

на

n = 3,

проверяем

знак

значения

f (xn

0,001)

f ( 10,260).

Так как

f ( 10,260) 0,

из этих чисел дает искомое приближение.

то

10,261 ξ 10,260,

и любое

Пример 3.5
Методом	золотого	сечения	найти	корень	уравнения
f (x) x 4 2x3 x 1 0 на отрезке [0;1] с точностью = 10-2.
Решение представлено в табл. 3.4.

Таблица 3.4

Реализация метода золотого сечения

шаг	ai	bi	f(ai)	f(bi)
0	0	1	-1	1
1	0.6180	1	-1	1
2	0.8541	1	-0,0758	1
3	0.8541	0,9443	-0,0758	0,5347
4	0.8541	0,9098	-0,0758	0,2817
5	0.8541	0,8885	-0,0758	0,1377
6	0.8541	0,8754	-0,0758	0,0534

	c	5 1	b a a
		2

0,6180

0,85414

0,9443

0,9098

0,8885

0,8754

0,8673

f(c)

-1,0001

-0,0758

0,5347

0,2817

0,1377

0,0534

0,003

Знак

f (c)

-

-

+

+

+

+

+

В качестве приближенного корня возьмем среднее арифметическое,

то есть

ab

x 6 6

2

0,8648

.

Пример 3.6
Графически решить уравнение	x lg(x) 1.

Решение. Запишем наше уравнение в виде равенства

lg(x)

1 x

. Отсюда

ясно, что корни исходного уравнения могут быть найдены как абсциссы

точек пересечения логарифмической кривой y lg(x) и гиперболы

y

1 x

.

Построив эти кривые, приближенно найдём один корень Графическое решение представлено на рис. 3.20.

x	0

2,5

.

91

Рис. 3.20. Графическое решение уравнения

Пример 3.7

С помощью графического метода отделить корни трансцендентного уравнения и уточнить их методом Ньютона с точностью =0,00001:

	ln(x) (x 1)	2	0,15	0 .
Решение.	Запишем наше уравнение в			виде	y1	y2	, где	y1	ln(x) ;
y2 (x 1)2	0,15 . Строим графики данных функций.

Из рис. 3.21 видно, что данное уравнение имеет два корня: первый корень принадлежит отрезку [0.1; 1], а второй [1.1; 2].

Уточним корни методом касательных. Для этого вычислим производные

f		1	2 x 1 ;
	(x)	x

4

2

y1(x)

f

	(x)		1
		2	2
	x

.

y2(x)

0

1

2

3

2

4

x

Рис. 3.21. Графический метод отделения корней трансцендентного уравнения

Итерационная формула метода Ньютона в данном случае имеет вид

xn xn 1	ln(x	n 1		) (x		n 1	1)2 0,15		,
				1	2 x		1

xn 1 n 1

где n = 1, 2, 3, … .

92

f (0,1) 2.96	;	f		102	;	f (0,1) f
			(0,1)				(0,1)

Результаты вычислений представим

0 , поэтому х0=0,1. в виде табл. 3.5.

Таблица 3.5

Реализация графического метода (первый корень)

k	xk	f(xk)	f’(xk)
0	0,1	-2,962585	11,8
1	0,351067	-1,31789	4,146330
2	0,668912	-0,36172	2,157139
3	0,836598	-0,05511	1,522122
4	0,872805	-0,00222	1,400120
5	0,874392	-0,000001	1,394868
6	0,874395	0,0000001	1,394858
7	0,874395

Аналогично получаем результаты для второго корня в табл. 3.6

								Таблица 3.6
			Реализация графического метода (второй корень)
k			0	1	2	3		4
xk			2	1,8954315	1,8856575	1,8855667		1,8855667
f(xk)			-0,1568528	-0,0123510	-0,0001089	0,00000001
f x	k		-1,5	-1,2692786	-1,2409892	-1,2407889

3.8.Пример применения нелинейных уравнений

вэкономике

Одной из распространенных экономических задач является задача максимизации прибыли предприятия. Известно, что балансовая прибыль есть разница между выручкой и затратами на производство продукции P=N-ZВ. В общем случае выручка от реализации продукции может быть представлена полиномом второй степени от количества продукции N=a0Q+a1Q2. Нелинейность может быть связана с тем, что в условиях монополии цена единицы продукции k может уменьшаться с ростом количества выпущенной продукции Q:

k=a0+a1Q (a0>0,a1<0).

В свою очередь, функция затрат может быть представлена полиномом 3-й степени Z=b0+b1Q+b2Q2+b3Q3. Кубическая нелинейность может объясняться тем, что при производстве малой партии товаров издержки быстро растут, затем с ростом Q темп роста издержек уменьшается, но по достижении некоторого критического значения Q начинает работать «закон убывающей отдачи», в соответствии с которым

93

издержки вновь начинают расти ускоренными темпами. Прибыль максимальна, когда dP/dQ= 0.

Пример 3.8

С помощью пакета Excel решим задачу нахождения максимальной прибыли, полагая заданными коэффициенты: b0 = 10, b1=1, b2=-0,1, b3=0,01, a0=5, a1= -0,1, полагая также, что выпущен 21 вид продукции. Последовательность действий при реализации в пакете Excel (рис. 3.22).

Рис. 3.22. Реализация максимизации прибыли на предприятии в Microsoft Ecxel

1.Оформить заголовок в строке 1 «Максимизация прибыли».

2.В ячейки A3, ВЗ, СЗ, D3 и ЕЗ записать заголовки рядов соответственно Q, N, Z, P и dP/dQ.

3.В ячейки F3, F4, F5, F6, F9, F10 записать названия коэффициентов соответственно b0, b1, b2, b3, a0, a1.

4.В ячейки G3, G4, G5, G6, G9, G10 записать значения коэффициентов

соответственно 10; 1; -0,1; 0,01; 5; -0,1.

5.В ячейку Н5 ввести текст «Издержки Z=b0+bl Q+b2 Q^2+b3 Q^3».

94

6.В ячейку Н6 ввести текст «Выручка N=a0 Q+a1 Q^2».

7.В ячейку Н7 ввести текст «Прибыль P=N-Z».

8.В ячейки А4 и А5 ввести первые два значения аргумента 0 и 1.

9.Выделить ячейки А4:А5 и протащить ряд данных до конечного значения «21», убедившись в правильном выстраивании арифметической прогрессии.

10.В ячейку В4 ввести формулу «=A4*$G$9+A4*A4*$G$10».

11.Скопировать формулу на остальные элементы ряда, используя прием протаскивания. В интервале В4:В25 получен ряд результатов вычисления выручки N(Q).

12.В ячейку С4 ввести формулу

«=$G$3+A4*$G$4+A4*A4*$G$5+A4*A4*A4* $G$6».

13.Скопировать формулу на остальные элементы ряда, используя прием протаскивания. В интервале С4:С25 получен ряд результатов вычисления издержек Z(Q).

14.В ячейку D4 ввести формулу «=B4-C4».

15.Скопировать формулу на остальные элементы ряда, используя прием протаскивания. В интервале D4:D25 получен ряд результатов вычисления прибыли P(Q).

16.В ячейку Е4 ввести формулу «=($G$9-$G$4)+2*($G$10-$G$5)*A4-

3*$G$6* А4*А4».

17.Скопировать формулу на остальные элементы ряда, используя прием протаскивания. В интервале Е4:Е25 получен ряд результатов вычисления dP/dQ для различных значений Q.

18.Построить на одной диаграмме графики зависимостей N(Q), Z(Q) и P(Q), используя соответствующие ряды данных.

19.Построить на отдельной диаграмме зависимость dP/dQ от Q. Точка пересечения графика с осью абсцисс дает значение Q, соответствующее максимальной прибыли (шаговый метод).

95

Контрольные вопросы

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

2.Каким условиям должна соответствовать функция f(x) и что они гарантируют?

3.Для чего необходимо отделять корни?

4.В чем состоит метод половинного деления?

5.Выведите формулу погрешности для метода половинного деления.

6.Запишите расчетную формулу метода Ньютона и дайте геометрическую интерпретацию метода.

7.Как выбрать начальное приближение для метода Ньютона?

8.В чем сущность метода итерации?

9.Сформулируйте достаточные условия сходимости метода итерации?

10.Запишите расчетную формулу метода золотого сечения и дайте геометрическую интерпретацию метода.

96

4. Аппроксимация и интерполяция функций

			В	вычислительной			математике
		существенную			роль	играет	интерполяция
		функций (термин интерполирование впервые
		употребил Дж. Валлис в 1656 г. при составлении
		астрономических				и
		математических			таблиц),
		то	есть	построение		по
		заданной функции другой
		(как	правило,		более
		простой),			значения
	Дж. Валлис	которой		совпадают		со
	(1616 – 1703)
		значениями			заданной
	функции в некотором числе точек. Причем
	интерполяция имеет как		практическое, так			и
	теоретическое значение.					Жозеф Луи Лагранж
							(1736 – 1813)
	На практике часто		возникает задача			о

восстановлении непрерывной функции по ее табличным значениям, например, полученным в ходе некоторого эксперимента. Для вычисления многих функций оказывается эффективным приблизить их полиномами или дробно-рациональными функциями, для этого используют интерполяционный многочлен Лагранжа, который был предложен в 1797 г. в связи с решением задач дифференциального исчисления; интерполяционную формулу Ньютона, которая была впервые опубликована в «Методе разностей» (1711 г.). Также применяют кубический сплайн – это функция, являющаяся сравнительно новым математическим изобретением. Сплайн моделирует старое механическое устройство. Чертёжники издавна использовали длинные тонкие рейки из какого-нибудь упругого материала в качестве лекал, проводя с их помощью плавные кривые через заданные точки. Эти рейки или механические сплайны закрепляют, подвешивая грузила в точки интерполяции, называемые исторически узлами.

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

97

4.1. Конечные разности различных порядков
Пусть y=f(x) заданная функция. Обозначим через	x h

фиксированную величину приращения аргумента (шаг). Тогда выражение

y f ( x ) f ( x x ) f ( x )

(4.1)

называется первой конечной разностью функции у. Аналогично

определяются конечные разности высших порядков	∆ny =∆(∆n -1y)
(n= 2, 3, . . . ) .
Справедливо утверждение: если Pn(x) = a0xn + a1xn-1 + ... + an полином
n-й степени, то ∆nРn(x) = n!a0hn = const, где ∆x = h.
Действительно, имеем
∆ Pn(x) = Pn(x+h) - Pn(x) = a0[(x+h)n- xn] + a1[(x+h)n-1- xn-1]+ ... +
+ an-1[(x + h) – h].
Раскрыв по биному Ньютона круглые скобки, легко убедиться, что
∆Pn(x) представляет собой полином (п-1)-й степени:
∆Pn(x) = b0xn-1 + b1xn-2 + ... + bn-1,
где b0= nha0.
Рассуждая аналогично, приходим к выводу, что	вторая разность

∆2Pn(x) есть полином (п - 2)-й степени:

∆2Pn(x) = c0xn-2 + c1xn-3 + ... + cn-2,

причем с0 = (n - 1)hb0 = n(n - 1)h2a0.

Проводя последовательно аналогичные рассуждения, мы, в конце концов, установим, что ∆nPn(х) = n!a0hn= const. Как следствие получаем ∆s Рп (х) = 0 при s >п. Символ ∆ (дельта) можно рассматривать как оператор, ставящий в соответствие функции y = f(x) функцию ∆y = f(x + ∆x) – f(x) (∆x постоянно). Легко проверить основные свойства оператора ∆:

1)∆(u + v) = ∆u +∆v,

2)∆(Cu) = C∆u (С постоянная),

3)∆m(∆ny) = ∆n+my,

где т и п целые неотрицательные числа, причем по определению полагают ∆nу = у.

Из формулы (4.1) имеем f(x + ∆x) = ∆f(x) + f(x).

Отсюда, рассматривая ∆ как символический множитель, получим

f(x + ∆x) = (1 +∆)f(x).	(4.2)
Последовательно применяя это соотношение п раз, будем иметь
f(x +n ∆x) = (1 +∆)nf(x).	(4.3)

Воспользовавшись формулой бинома Ньютона, окончательно выводим

98

где

C	m		n(n 1)...[n (m 1)]
	n		m!

n		f (x) ,
f (x n x) Cn
m	m
m 0

число сочетаний из n элементов по m.

(4.4)

Таким образом, с помощью формулы значения функции f(х) выражаются через различных порядков.

Воспользовавшись тождеством

(1 ) 1

и применяя бином Ньютона, получаем

(4.4) последовательные ее конечные разности

(4.5)

n

f (x) [(1 ) 1]

n

f (x) (1

)

n

f (x) C

1

(1

)

n 1

f (x) C

2

(1

)

n 2

f (x)

n

n

... ( 1)

n

f (x).

Отсюда в силу формулы (4.3) будем иметь

	f (x) f (x n x) C		f [x (n 1) x] C		f [x (n 2) x] ... ( 1)		f (x).
n		1		2		n
		n		n

(4.6)

Формула (4.6) дает выражение конечной разности п-го порядка функции f(х) через последовательные значения этой функции.

Пусть функция f(х) имеет непрерывную производную f(n)(x) на отрезке (х , x+n∆x). Тогда справедлива формула

∆nf(x) = (∆x)nf(n)(x + θn∆x),

(4.7)

где 0<θ<1.

Формулу (4.7) проще всего доказать, используя метод математической индукции. В самом деле, при п = 1 мы получаем теорему Лагранжа о конечном приращении функции и, следовательно, формула

(4.7) верна. Пусть теперь при k < п имеем

∆kf(x) = (∆x)kf(k)(x + θ’ n∆x),

где 0<θ’<1.

Тогда	f (x)	[ f (x x) f (x)] (x)		[ f		(x x θ k x) f		(x θ k x)].
k 1	k		k		(k )		(k )

Применяя теорему Лагранжа к получившемуся приращению производной f(k)(x), будем иметь ∆k+f(x) = (∆x)k+1f(k+1)(x + θ’k∆x+ θ”∆x), где 0 < θ”< 1. Полагая

	θ k θ	θ ,	(4.8)
	k 1

окончательно получим ∆k+1f(x) = (∆x)k+1f(k+1)(x+θ(k+1)∆x), причем, очевидно, 0 <θ< 1.

Таким образом, установлен переход от k к k+1 и, следовательно, формула (4.7) доказана.

Из формулы (4.7) имеем

f (n) (x θn x) n f (x) . ( x)n

99

Отсюда, переходя к пределу при ∆x→0 и предполагая, что производная f(n)(х) непрерывна, получим

			n	f (x)
f	(n)	(x) lim			.
			(x)
		x 0
				n

(4.9)

Следовательно, при малых ∆x справедлива приближенная формула

			n	f (x)
f	(n)	(x)			.
			(x)
				n

(4.10)

4.2. Таблица разностей

Часто приходится рассматривать функции y = f(x), заданные табличными значениями yi=f(xi) для системы равноотстоящих точек xi

(i = 0, 1,2, . . . ) , где ∆ xi = xi+1 - xi =h= const.

Конечные разности последовательности yi естественно определяются соотношениями

∆yi= yi+1 - yi,

∆2yi=∆(∆yi) = ∆yi+1 - ∆yi ,

. . .

∆nyi=∆(∆n-1yi) = ∆n-1yi+1 - ∆n-1yi.

Из первого равенства имеем yi+1 = yi +∆yi=(1+∆)yi. Отсюда последовательно выводим

yi+2 = (1+∆) yi+1 =(1+∆)2 yi , yi+3 = (1+∆)yi+2 =(1+∆)3 yi ,

. . .

yi+n = (1+∆)n yi .

Использовав формулу бинома Ньютона, получим

y	y	1	y	2	2	n
		C		C	y	... y .
i n	i	n	i	n	i	i

Обратно, имеем

	n y [(1 ) 1]n y		(1 )n y	C1	(1 )n 1 y			C2 (1 )n 2 y		... ( 1)n y ,
	i	i	i	n			i	n	i	i
или	n	1	2	...	( 1)	n	yi .
	yi	yn i Cn yn i 1 Cn yn i 2
	Например, ∆2yi= yi +2 - 2yi+1+yi, ∆3yi= yi +3 - 3yi+2+3yi+1 - yi и т.д.
	Заметим, что для вычисления n-й разности ∆nyi								нужно знать п + 1
членов yi, yi+1, ... , yi+n данной последовательности.

100