Практикум j практическое примщенше численных методов

xRootl = bisection(f, 2, 4) print ’root(l) = xRootl xRoot2 = bisection(f, 6, 8) print ’root(2) = ’, xRoot2 xRoot3 = bisection(f, 8, 10)


Рис. 7.1 График функции у = (1 + х²)е ^х + sin(:r) print ’root(3) = xRoot3



На интервале [0,10] найдены три корня (см. рис. 7.1).
Упражнение 7.2 Напишите программу для нахооюдения решения сисгпе- мы нелинейных уравнений F(x) = 0 методом Ньютона при численном вы- числении матрицы Якоби. С ее помощью найдите приблиэюениое решение системы
(3 + 2xi)xi - 2х₂ = 3,
(3 + 2Xi)xi - Xi_i - 2x_i+1 =2, i = 2,3,..., n - 1,
(3 + 2z_n)z_n -x_n-i = 4
при n = 10 и сравните его с точным решением ж» = 1, г = 1,2,...,п.
В модуле newton функция jacobianO предназначена для вычисления эл^с' ментов матрицы Якоби с использованием разностной производной. Функи!¹*¹ newton () реализует метод Ньютона при решении соответствующей систем^ линейных уравнений с помощью LU-разложения (функция solveLUO из мо­дуля lu).
Л■ bit::;: i::;: ■::: I- -
import numpy as np import math as mt from lu import solveLU def jacobian(f, x):
и и и
Calculation of the Jacobian using finite differences.
к и и
h = 1.0e-4 n = len(x)
Jac = np.zeros((n,n), ’float’) fO = f(x)
for i in range(n): tt = x[i] x[i] = tt + h fl = f(x) x[i] = tt
Jac [:,i] = (fl - f0)/h return Jac, fO
def newton(f, x, tol=1.0e-9):
и и n
Solves the system equations f(x) =0 by
the Newton method using {x} as the initial guess.
Solve the linear system Ax = b by lu module.
и и и
iterMax = 50
for i in range(iterMax):
Jac, fO = jacobian(f, x) if mt.sqrt(np.dot(fO, fO) / len(x)) < tol: return x, i dx = solveLU(Jac, fO) x = x - dx
print ’Too many iterations for the Newton method’
В тестовой задаче метод Ньютона применяется при начальном приближении х9 = 0, г = 1,2, ...,п и с использованием следующей программы.
import numpy as np from newton import newton n = 10 def f(x):
f = np.zeros((n), ’float’) for i in ranged,n-1):
f[i] = (3 + 2*x[i])*x[i] - x[i-1] - 2*x[i+l] - 2 f CO] = (3 + 2*x[0] )*x[0] - 2*x[l] - 3 f[n-l] = (3 + 2*x[n-l])*x[n-l] - x[n-2] - 4 return f
xO = np.zeros((n), ’float’) x, iter = newton(f, xO) print ’Newton iteration = ’, iter print ’Solution:\n’, x
1ИШв!ШИ11Н
111
При используемых параметрах метода точное решение задачи совпадает с численным.

4. 
   Задачи
   

Задача 7.1 Напишите программу для решения нелинейного уравнения f(x) = 0 методом секущих. Используйте ее для решения уравнения 4 sin(a;) + 1 — х = 0 на интервале /-10,10/.
Задача 7.2 Напишите программу для решения нелинейного уравнения f(x) = 0 методом Ньютона. С ее помощью найдите полоэюительные корни уравнения х² — lOsin(rc) = 0.
Задача 7.3 Напишите программу для нахоэюдения решения системы нели­нейных уравнений F(x) = 0 методом Ньютона при аналитическом задании матрицы Якоби. С использованием этой программы найдите приблиэюеинос решение системы уравнений
п
n-£ (cos Xj + г(1 - cos rr*) - sin я*) = 0 г = 1,2,n,
j=i
при n = 20 и начальном приблиэюепии х® = 1/п, г = 1,2,...,п.
Задача 7.4 Напишите программу для нахоэюдепия решения системы нели­нейных уравнений F(x) = 0 методом Зейделя при решении нелинейного урав­нения методом Ньютона. Рассмотрите возмоэюности использования этой программы при приблиэюенпом решении системы уравнений из упраэюпс- ния 7.2.
Задачи минимизации функций
Среди основных проблем вычислительной математики можно отметить зада­чи минимизации функций многих переменных (задачи оптимизации). Поиск минимума часто проводится при некоторых дополнительных ограничениях — условная оптимизация. Для численного решения таких задач используют­ся итерационные методы. В задачах с ограничениями применяются методы штрафных функций. Простейшей задачей рассматриваемого класса является поиск минимума одномерной функции.

Основные обозначения
X {Я'г} {*^1? ? • • • > Хп }	— n-мерный вектор
/0*0	— функция одной или п переменных
f(x) —> min, х Е R^n	— задача безусловной оптимизации
f(x) —> min, x E X	— задача условной оптимизации
X	— допустимое множество
x^k n	— приближенное решение на к-й итерации
(x,y)=Y,^X*yi	— скалярное произведение
1=1

1. 
   Поиск минимума функции многих переменных
   

Для заданной функции /(ж), определенной на допустимом множестве X из евклидова пространства Rⁿ, ищется точки минимума (максимума) функции /(ж), т.е.
f(x) —> min, х Е X. (8.1)
^очка х* Е X есть точка глобального минимума функции f(x) на множестве если
/(я*) < /(ж), Vt Е X, (8.2)
^и точка локального минимума, если f(x*) < f(x) в окрестности точки х* Е X.
Задача (8.1) называется задачей безусловной оптимизации, если X = Rⁿ, т.₀
f{x) —> min, х е Rⁿ. (8.3)
Если X некоторое подмножество пространства Rⁿ, то мы имеем задачу услоц, ной оптимизации. Такие задачи существенно сложнее для численного реше­ния, чем задачи безусловной минимизации. Ограничения могут формулиро- ваться в виде равенств (например, ffi(x) = 0, г 1,2,... , га) или неравенств (9i(x) <0, г = 1,2,... ,тп).

2. 
   Методы решения задач оптимизации
   

Вычислительные алгоритмы для приближенного решения задачи оптими­зации чаще всего строятся на использовании необходимых и достаточных условиях оптимальности, т.е. на условиях, которые имеют место в точке ми­нимума. Реализация такого подхода связана с решением соответствующих нелинейных уравнений итерационными методами.
Поиск минимума функции одной переменной
Пусть X = [а,Ь] и кусочно-непрерывная функция f(x) имеет в некоторой точке х* Е X один минимум. Мы отметим прежде всего простейшие итераци­онные методы решения, задачи минимизации, наиболее полно учитывающие специфику одномерных задач.
Вычислим функцию на концах отрезка и в двух внутренних точках х¹ и х² < х¹. Будем считать, что эти точки симметричны относительно середины отрезка [а, Ь]. В методе золотого сечения точки х¹ и х² выбираются так, чтобы отношение длины всего отрезка [а, Ь] к длине большей из его частей [а, х¹] равно отношению длины большей части [а, .т¹] к длине меньшей части [х¹^]:
b — а х¹ — а х¹ — а х² — а
Далее проводится сравнение значений функции в четырех точках а,ж²,х^д,^ и выберем точку, в которой значение функции наименьшее. Пусть это будет точка х², тогда минимум функции достигается в одном из прилегающих к этой точке отрезков: [а, ж²] или х²,х¹] и поэтому в дальнейшем можно рас­сматривать проблему минимизации на отрезке [а,я¹]. После этого процесс повторяется — в соответствии с правилом золотого сечения делится точкой х³ отрезок [а,ж¹].
Для минимизации функции одной переменной широко используются методы полиномиальной интерполяции. В этом случае с использованием ранее най­денных точек строится интерполяционный полином, точка минимума кото­рого принимается за очередное приближение. В методе парабол используется интерполяционный многочлен второго порядка.
Пусть, например, известны три приближения х^к~² < х^к~¹ и х^к~² < х^к < х^к~¹, причем f(x^k) < f(x^k~²) и f(x^k) < f(x^k~^l). Новое приближение ищется как решение задачи минимизации
Ь₂(х) —> min, (8.5)
где L₂ (х) — интерполяционный многочлен второго порядка, построенный по узлам х^к~²,х^к~¹ и х^к.
Интерполяционный многочлен Ньютона имеет вид
L₂(x) = f{x^k) + {х - а;*^:)/(ж'^г,а;'^г_1)+
+
где


f(x^k,x^k *) =


f(x^k) - /(х^¹)