Численные методы

Преобразуем уравнение к виду: y f (x, y) xy x12 . Найдем числа:

k1	h f (x, y) 0,2		0		1			0,2;
			1			2
					1

k	h f (x h ,y k1) 0,2 0,1												1						0,18;
													2
2	2	2
								1,1			1,1
k	h f (x h ,y		k2		) 0,2 0,09								1					0,18;
3	2	2												2
							1,1						1,1
k	h f (x h,y k ) 0,2							0,18					1						0,17.
													2
4		3						1,2
											1,2
Следовательно,		у		0		1 k		2k	2		2k						3		k	4	0,18,	,
						6	1

т. е. y1 y0 y0 0 0,18 0,18

Аналогичным образом находим

	k h f (x,y) 0,2 0,18					1		0,17;
					2
	1		1,2
					1,2

k	h f (x h	,y k1) 0,2 0,26										1				0,15;
												2
2	2	2					1,3
									1,3
k	h f (x h	,y	k2	) 0,2 0,25								1				0,15;
												2
3	2	2					1,3
											1,3
k	h f (x h,y k ) 0,2 0,33											1				0,14.
												2
4		3
						1,4			1,4
Следовательно, y1					1	(k1	2	k2				2 k3 k4 ) 0,15,
					6

т. е. y2 y1 y1 0,18 0,15 0,33 и т. д.

Рассмотрим решение данного примера в Mathcad:

50

Рис. 36. Пример реализации метода Симпсона в Mathcad

5.4. Решение дифференциальных уравнений с помощью встроенных функций Mathcad

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

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

Для решения задачи Коши MathCad предлагает несколько функций, из которых чаще всего используется rkfixed, реализующая метод Рунге-Кутта с фиксированным шагом.

51

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

Этапы решения дифференциального уравнения:

1. Задать начальное значение искомой функции в качестве переменной y0 (для уравнения первого порядка) или вектор (для уравнений более высокого порядка или систем дифференциальных уравнений).

2. Разрешить ОДУ относительно производной, т. е. привести к виду y f (x, y) . Полученную функцию f (x, y) вписать в документ в виде

D(x, y): f (x, y0 ).

Если уравнение второго порядка, то его нужно разрешить относительно старшей производной, применить подстановку y t и записать систему ОДУ первого порядка. Начальные условия, как и функция D(x, y), записываются в виде вектора.

3. Далее производится решение с помощью функции rkfixed, реализующей метод Рунге – Кутта четвертого порядка с фиксированным шагом.

Функция rkfixed(y, x1, x2, npoints, D) имеет следующие параметры :

y – вектор начальных значений размерности n (для уравнения первого порядка n 1);

x1,x2 – граничные точки интервала, на котором ищется решение; npoints – количество точек, в которых ищется решение;

D – вектор из функций правых частей уравнения.

	Ниже показано решения уравнения y	y x	с начальным условием
		y x
	у(0) = 1.

52

Рис. 37. Решение ОДУ в Mathcad с помощью функции rkfixed

Имеется также функция Rkadapt (у, x1, x2, npoints, D), которая проверяет насколько быстро меняется решение, и соответственно адаптирует величину шага интегрирования. Представлен и ряд других функций для решения так называемых «жестких систем». Результат применения этих функций – матрица, содержащая n + 1 столбец. Первый столбец Y<1> содержит координаты узлов сетки, i-ый столбец Y<1> – вычисленные приближенные значения решения уi-1(x), (i = 2, …, n + 1) в узлах сетки. Пример применения функции rkfixed приведен ниже.

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

В Mathcad имеется также удобная конструкция Given…оdesolve, предназначенная в основном для решения линейных уравнений. Внутри этой конструкции дифференциальные уравнения и начальные условия пишутся фактически в естественном виде. В самой функции оdesolve (х, b) в качестве параметров указана независимая переменная х и конечные значения b, вплоть до которого будет строиться решение.

53

5.5. Индивидуальные задания для лабораторной работы № 4

Найти решение дифференциального уравнение по алгоритму:

1.Приведете уравнение к нормальному виду, задайте исходные данные для решения задачи Коши, и вычислите шаг интегрирования.

2.Используя микрокалькулятор, найдите решение методами РунгеКутта и Эйлера.

3.Используя средства программирования системы Mathcad, составьте программу численного решения задачи Коши дифференциального уравнения методами Эйлера и Рунге-Кутта. Постройте графики.

4.Решите задачу Коши с использованием встроенной функции rkfixed. Постройте график.

5.Убедитесь, что это решения совпадают.

6) Сформулируйте выводы по проделанной работе.

Из таблицы выбрать свой вариант задания (y f (x,y), x [a,b], y0 y(a))

№			f(x,y)												a	b	y(a)	№	f(x,y)															a	b	y(a)
1	x sin(							y				)			1	3,4	2	16	y sin(						y			)						2	2,8	1
								x																	x
2		x e y													3	3,8	4	17	x y y															2	4,4	3
3		x				y									0	1,6	2	18	y		sin x													3	3,8	1
					ln y
																						y
4			1			y									0	3,2	1	19							y									1	2,6	2
																			x cos(											)
		x y
																									x
5	y		x cosy												1	1,8	1	20	x cos(					y					)					3	5,4	1
																									x
6	x			1					x						3	5,4	2	21		x y														2	3,6	2
									y											x y
7		e x y				y									2	5,2	1	22	y e x y															1	3,4	2
8			y x 2												1	1,8	4	23	x e y x															3	3,8	1
9		x		sin y											0	1,6	1	24	x cos(				x										)	2	5,2	1
						y																x y
10							y								4	6,4	0,5	25									x							3	4,6	1
		y																	1 cos(													)
																										y
							x
11	y ln(y x)														0	1,6	2	26	y 2 x 2															2	2,8	3

54

	12		y e y x				3	4,6	0	27	y		x						4	4,8	1
													ln(x)
	13		y 2		x		1	1,8	2	28	ln(		x				)		1	2,6	1
												x y
	14		ln x		y		4	5,6	3	29	x sin(			y		)			2	2,8	1
			x											x
	15	y			1		3	4,6	1	30	cos(x				y)				0,6	2,2	0
				x y									x

5.6.Вопросы для самоконтроля

1.Что называется обыкновенным дифференциальным уравнением? Записать в нормальной форме дифференциальное уравнение первого порядка.

2.Дайте определение частного и общего решения дифференциального уравнения.

3.Как ставится задача Коши для дифференциального уравнения первого порядка? При каких условиях эта задача имеет, причем единственное, решение?

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

5.В каком виде представляется решение задачи Коши, найденное численным методом?

6.Записать расчетную формулу метода Эйлера. Каков геометрический смысл метода Эйлера?

7.Записать расчетные формулы для метода Рунге – Кутта.

8. Записать формулы для расчета погрешности в методах Эйлера и Рунге – Кутта.

55

6.ЧИСЛЕННЫЕ МЕТОДЫ В ЭКОНОМИЧЕСКИХ ЗАДАЧАХ

6.1.Нелинейные уравнения в экономических задачах

Задача 1. Облигация сроком n лет с выплатой процентов раз в год в конце года по норме g куплена по курсу K. Определить доходность к погашению i. Для этого необходимо решить уравнение f(i) = 0 вида:

	(1 i) n g	1 (1 i) n			K	0
		i		100

(Например, n = 8; g = 0,12; K = 101).

Задача 2. Сумма наращенных инвестиций составляет S млн руб. при первоначальной цене P тыс. руб. Инвестиции окупятся за n лет. Определите процентную ставку i, по которой начисляются проценты на долг. Для этого необходимо решить уравнение f(i) = 0 вида:

PS i (1 i)n 1 0

(Например, n = 9; S = 1,5; P = 0,11).

Задача 3. Создаётся фонд, равный S млн руб. На протяжении n лет раз в конце каждого года производятся платежи в размере R тыс. руб. На собранные средства начисляются проценты. Необходимо найти годовую ставку процента i. Для этого необходимо решить уравнение f(i) = 0 вида:

S R (1 i)n 1 0 i

(Например, n = 5; S = 2,4; R = 0,15).

Задача 4. Кредит выдается при ставке i процентов годовых. Для его погашения создан фонд путем вложения равных сумм. В фонде начисляются g процентов годовых. Необходимо погасить основной долг и проценты. Определить годовой взнос ie, если в фонд вносятся деньги в начале года в течение n лет. Проценты погашаются в конце каждого года средствами, накопленными в фонде. Для этого необходимо решить уравнение f(ie) = 0 вида:

56

1 (1 i ) n

i ie (1 (1 e i) n ) (1 g) 0 (Например, n = 7; i = 0,19; g = 0,1).

6.2. Интерполирование функций в экономических задачах

Задача 1. Для значений f(ti) экономического показателя деятельности фирмы в моменты времени ti [t0 ,tn ], i 0...n , требуется вычислить значение экономического показателя в момент времени x с помощью многочлена Лагранжа. Определить вычислительную погрешность многочлена Лагранжа на интервале [t0 ,tn ], считая, что f(ti) заданы с верными знаками. Построить на одном графике значения многочлена Лагранжа и исходных данных, а на другом – вычислительную погрешность на интервале [t0 ,tn ]. Проанализировать полученные результаты и сделать выводы.

ti		1	2,2	3,4	5,4	6	6,7	8,1	x=3,8
f(ti)		-5,121	-1,346	2,75	5,011	5,771	8,321	12,841
	Задача 2. Для значений f(ti) экономического показателя деятельности
фирмы в моменты времени ti [t0 ,tn ],						i 0...n , требуется вычислить зна-

чение экономического показателя в момент времени x с помощью многочлена Ньютона. Определить вычислительную погрешность многочлена Ньютона на интервале [t0 ,tn ], считая, что f(ti) заданы с верными знаками. Построить на одном графике значения многочлена Ньютона и исходных данных, а на другом – вычислительную погрешность на интервале [t0 ,tn ]. Проанализировать полученные результаты и сделать выводы.

ti	0,7	1,2	3	3,9	4,2	5,8	7,7	x=0,9
f(ti)	-20,242	-15,302	-2,789	0,356	4,258	10,326	13,249

57

6.3. Определенные интегралы в экономических задачах

Пусть функция интенсивности спроса имеет вид s(t) k e c z(t),

где k – коэффициент, задающий интенсивность спроса; c – цена единицы товара (руб), z(t) –объем товара на рынке (руб.). Объем востребованного товара за время T равен значению интеграла

I T s(t)dt .

0

Кроме интенсивности спроса на востребованность товара, существенное влияние оказывает реклама. Будем предполагать, что влияние рекламы описывается убывающей функцией

f (t,x)	R		0,5x2 5,75x T 2	,
	T		(t 1) (x 1)

линейно зависящей от объема средств R (руб.), затрачиваемых на рекламу. Таким образом, в начальный момент времени реклама оказывает наибольшее воздействие, которое постепенно убывает. При этом объем востребованного товара равен значению интеграла

T (t )

I k e c z(t) f (t,x)dxdt

00

Задача 1. Определить общее количество востребованного товара в зависимости от интенсивности спроса s(t) на временном интервале [0, T] по формуле левых прямоугольников. Для вычисления интеграла с точностью ε 0,5 10 5 использовать правило Рунге. Вывести значение интеграла и трудоемкость метода, которая определяется числом вычислений подынтегральной функции. При

T 5; k 12; c 1,5; z(t) sin(Tt ); ψ(t) (t T1)3 ; R 20

58

Задача 2. Определить общее количество востребованного товара в зависимости от интенсивности спроса s(t) на временном интервале [0, T] по формуле Симпсона. Для вычисления интеграла с точностью ε 0,5 10 5 использовать правило Рунге. Вывести значение интеграла и трудоемкость метода, которая определяется числом вычислений подынтегральной функции при следующих значениях:

T 8; k 11; c 3,5;

z(t)

ln(t 1); (t) (t 1) e t ;

R 21

6.4.Обыкновенные дифференциальные уравнения

вэкономических задачах

Задача 1. Решить дифференциальное уравнение, описывающей динамику односекторной экономики на временном отрезке [Т,0], где Т = 2 года, с шагом h, равным одному месяцу методом Эйлера. Для моделирования использовать производственную функцию Кобба – Дугласа с параметрами: A 1; K (0) 10; L 10, α – заданное значение, β = 1 – α , Параметр µ задать таким образом, чтобы амортизация основных фондов к концу времени моделирования была равна 20 %. Оценить погрешность решения с помощью оценки Рунге. Моделирование осуществить для двух значений пара-

метра s:s1 (0;	0,5) , когда норма накопления меньше нормы потребления
и s2 (0,5; 1),	когда норма накопления больше нормы потребления. Выве-

сти графические изображения функций K(t), I(t), C(t) для значений s1 и s2. Проанализировать полученные результаты и сделать выводы при следую-

щих значениях: α = 0,25; s1 = 0,45; s2 = 0,8.

59