Нелиненйые уравнения

Цель лекции: Рассмотреть базовые методы нахождения приближенных решений нелинейных трансцендентных уравнений. Реализовать на C# основной метод.

Решение нелинейных уравнений является одной из самых распространенных задач математики. В лекции, посвященной особенностям вычислительных процедур, мы рассматривали простейший метод решения скалярного уравнения для вычисления корня квадратного из $$2$$. Этот метод - половинного деления имеет много достоинств. Главным достоинством этого метода является то, что для приближенных решений автоматически получается оценка точности. Второе достоинство данного метода состоит в том, что этот метод легко реализуем. Однако у этого метода есть и недостатки. Для реализации этого метода мы должны знать отрезок, на концах которого функция принимает значения разного знака. А также этот метод является сугубо скалярным и не может быть обобщен на многомерный случай.

Сначала рассмотрим задачу о нахождении корня одного уравнения$$f(x)=0.$$ Если предположить, что функция $$f(x)$$ является дифференцируемой, и ее производная есть непрерывная функция $$f'(x)$$, то находить приближенные решения возможно с помощью метода Ньютона. Формулы метода Ньютона получаются из следующих соображений. Пусть $$x^*$$ - корень функции $$f$$, а $$x_n$$ некоторое приближение к этому корню. Тогда имеем$$f(x^*)=f(x_n+(x^*-x_n))=f(x_n)+(x^*-x_n)f'(\xi)=0.$$ Если вместо неизвестной точки $$\xi$$ взять значение $$x_n$$, то следующее приближение можно вычислять по формулам$$x_{n+1}=x_n-\frac{f(x_n)}{f'(x_n)}.$$ Можно показать, что если начальное приближение $$x_0$$ $ достаточно близко к искомому корню, то итерации, вычисленные по методу Ньютона сходятся к этому корню. К сожалению, это условие носит скорее теоретический характер, поскольку не понятно, что значит "достоаточно близко" и как это условие проверить. Другим достаточным условием сходимости метода Ньютона является условие$$|f(x)f'(x)|<(f'(x))^2.$$ Также из геометрических соображений можно получить, что итерации сходятся к корню с той стороны, с которой выполнено условие$$f(x)f''(x)\ge0,$$ если функция $$f$$ дважды дифференцируема.

Метод Ньютона сходится достаточно быстро, однако у него есть недостатки, состоящие в том, что факт сходимости метода и его скорость сходимости зависят от начального приближения. Вторым серьезным недостатком этого метода является то, что для его реализации нужно знать аналитическое значение производной.

Для устранения последнего недостатка был разработан метод секущих, который иногда называется методом хорд. В этом методе производная заменяется конечной разностью$$f'(x_n)\approx\frac{f(x_n)-f(x_{n-1})}{x_n-x_{n-1}}.$$ Тогда итерации по методу секущих вычисляются по формулам$$x_{n+1}=x_n-\frac{(x_n-x_{n-1})f(x_n)}{f(x_n)-f(x_{n-1})}.$$ Для реализации метода секущих необходимо знать два предыдущих значения $$x_n$$ и $$x_{n-1}$$. Такие методы называются двухшаговыми. Метод секущих сходиться медленнее, чем метод Ньютона, зато этот метод не требуем знания производной функции $$f(x)$$.

Реализуем класс для вычислений по методу секущих.$$\begin{verbatim} class TSolve { int Max_Iter; double Eps; TRFunc f; public TSolve(int Max_Iter, double Eps, TRFunc f) { this.Max_Iter = Max_Iter; this.Eps = Eps; this.f = f; } public double Solve(double x0, double x1) { double xn1 = 0; double xn = x1; double xn_1 = x0; for (int n = 2; n <= Max_Iter; n++) { if (Math.Abs(f.CalcY(xn) - f.CalcY(xn_1)) < Eps) { break; } xn1 = xn - ((xn - xn_1) * f.CalcY(xn)) / (f.CalcY(xn) - f.CalcY(xn_1)); Console.WriteLine("x_{0} = {1}", n, xn1); xn_1 = xn; xn = xn1; } return xn1; } } \end{verbatim}$$ В этом классе мы используем класс $$TRFunc$$, реализующий числовые функции. Для начала вычислений необходимо задать два значения $$x_0$$ и $$x_1$$. При инициализации класса следует определить два параметра - максимальное количество итераций $$Max\_Iter$$ и минимальную невязку $$Eps$$. Вычисления останавливаются в двух условиях - либо когда превышено количество итераций либо, когда выполнится условие$$|f(x_n)-f(x_{n-1})|<Eps.$$ Проведем тестирование нашего метода.$$\begin{verbatim} TSinFunc Sin = new TSinFunc(-100, 100); TSolve Solve = new TSolve(10, 1E-12, Sin); Console.WriteLine("x_0 = {0}\tx_1 = {1}", 1, 1.5); Solve.Solve(1, 1.5); Console.WriteLine("\n\nx_0 = {0}\tx_1 = {1}", 1.7, 1.8); Solve.Solve(2, 2.2); \end{verbatim}$$ Мы будем вычислять корни функции $$f(x)=\sin x$$. В результате двух запусков нашего класса мы получим следующие результаты$$\begin{verbatim} x_0 = 1 x_1 = 1.5 x_2 = -1.69660749346521 x_3 = -0.102640761680463 x_4 = 0.0809387254524748 x_5 = -3.00669920823837E-05 x_6 = 3.28414216333597E-08 x_7 = -4.94283429392287E-18 x_8 = 7.70371977754894E-34 x_0 = 1.7 x_1 = 1.8 x_2 = 3.80414325114683 x_3 = 3.111014449759 x_4 = 3.14383344379986 x_5 = 3.14159232994785 x_6 = 3.14159265359006 x_7 = 3.14159265358979 \end{verbatim}$$ Видно, что если рассматриваемая функция имеет несколько корней, то различные значения начальных приближений могут приводить к различным решениям нашего уравнения.

Рассмотрим теперь вопрос о решении системы трансцендентных уравнений. Пусть необходимо решить следующую систему уравнений$$f_k(x_1,x_2,\dots,x_n)=0,\quad 1\le k\le n.$$ Как правило такие системы решают методом простых итераций. Для этого необходимо систему 13.1 записать в виде$$x_k=\varphi_k(x_1,x_2,\dots,x_n),\quad 1\le k\le n.$$ Или в векторной форме$$x=\varphi(x),$$ где $\varphi:\Bbb{R}^n\to\Bbb{R}^n$. Вместо нахождение корня системы 13.1 мы пришли к задаче нахождения неподвижной точки.

Метод простых итераций выглядит довольно просто. Необходимо задать начальное приближение $$x^0$$, а дальнейшие вычисления проводить по схеме$$x^{n+1}=\varphi(x^n).$$ Однако успех этого метода зависит от свойств функции $$\varphi$$ и начального приближения. Приведем одно достаточное условие сходимости метода простых итераций для случая, когда функции $$\varphi_k$$ являются непрерывно дифференцируемыми. Построим матрицу $$M$$ с элементами $$m_{ij}$$, определяемыми по формулам$$m_{ij}=\max\left|\frac{\partial \varphi_i}{\partial x_j}\right|.$$ Достаточным условием сжимаемости отображения $\varphi$ является условие$$\|M\|<1$$ в какой-либо матричной норме. В различных (но эквивалентных) нормах это условие выглядит так$$\sum\limits_{j=1}^nm_{ij}<1$$ или$$\sum\limits_{i=1}^nm_{ij}<1$$ или $$\sum\limits_{i,j=1}^nm^2_{ij}<1$$.

Ключевые термины

Корень функции - точка, в которой функция принимает нулевое значение.

Метод Ньютона - основной итерационный метод для нахождения приближенного решения системы нелинейных уравнений.

Начальное приближение - начальное значение последовательности, которая рассчитывается с помощью рекуррентной процедуры.

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

Краткие итоги: Рассмотрены численные методы нахождения приближенных корней трансцендентных уравнений. Приведена объектно-ориентированная реализация наиболее популярных методов нахождения корней уравнений.

Практическое занятие "Решение линейных и нелинейных уравнений"

Цель занятия

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

Практическая задача

Как мы уже отмечали на лекциях, важнейшим понятием для решения системы линейных алгебраических уравнений является число обусловленности матрицы. Вычислим это число для некоторых матриц. Согласно определению находить это число затруднительно. Однако для симметричных положительно определенных матриц это число может быть найдено с помощью собственных значений этой матрицы. Пусть симметричная положительно определенная матрица $$A$$ имеет собственные значения $$\lambda_{min},\ldots,\lambda_{max}$$, где $$\lambda_{min}$$ - минимальное собственное значение, а $$\lambda_{max}$$ - максимальное собственное значение. Тогда число обусловленности этой матрицы можно вычислить по формуле$$\mu(A)=\frac{\lambda_{max}}{\lambda_{min}}.$$

Рассмотрим матрицу$$A=\left(% \begin{array}{cc} 1 1/2 \\ 1/2 2 \\ \end{array}% \right)$$ Эта матрица является симметричной положительно определенной матрицей. Вековое уравнение для этой матрицы имеет вид$$\det(A-\lambda I)=\lambda^2-3\lambda+\frac{7}{4}=0.$$ Вычисляя собственные функции этой матрицы, мы имеем $$\lambda_{min}=\frac{3-\sqrt{2}}{2}$$, $$\lambda_{max}=\frac{3+\sqrt{2}}{2}$$. Число обусловленности этой матрицы равно$$\mu(A)=\frac{3+\sqrt{2}}{3-\sqrt{2}}=2.7836114...$$

Рассмотрим применение метода Холецкого для этой матрицы. Напомним, что метод Холецкого позволяет найти такую матрицу $$L$$, которая имеет вид$$L=\left(% \begin{array}{cc} l_{11} 0 \\ l_{21} l_{22} \\ \end{array} \right)$$ и для которой выполнено равенство$$A=LL^T.$$ Согласно методу Холецкого коэффициенты $$l_{ij}$$ находятся следующим образом$$l_{11}=\sqrt{a_{11}}=1,$$ $$l_{21}=\frac{a_{21}}{l_{11}}=\frac{1}{2},$$ $$l_{22}=\sqrt{a_{22}-l_{21}^2}=\sqrt{2-\frac{1}{4}}=\frac{\sqrt{7}}{2}.$$ Таким образом мы имеем матрицу$$L=\left(% \begin{array}{cc} 1 0 \\ 1/2 \sqrt{7}/2 \\ \end{array}% \right)$$

Теперь рассмотрим практические вопросы реализации метода Ньютона для нахождения решений трансцендентных уравнений. Напомним, что метод Ньютона применим для решения уравнений вида$$f(x)=0,$$ где $$f(x)$$ - непрерывно дифференцируемая функция. Основной вопрос, который здесь возникает, состоит в том, чтобы быть уверенным в сходимости итераций Ньютона. Основным достаточным условием сходимости итераций Ньютона является выполнение неравенства$$|f(x)f'(x)|<(f'(x))^2.$$ Найдем с помощью этого соотношения отрезок, на котором гарантированно сходится метод Ньютона для уравнения$$x^p=C,\quad p>0,\ C>0.$$ В этом случае мы имеем функцию $$f(x)=x^p-C$$. Производная этой функции имеет вид $$f'(x)=px^{p-1}$$. Тогда имеем неравенство$$(x^p-C)px^{p-1}<p^2x^{2p-2}.$$ Отсюда получаем$$x^p-C<px^{p-1}$$ $$-C<x^p(p-x).$$ Получаем, что на интервале $$(0,p)$$ метод Ньютона сходится.

Лабораторная работа "Решение линейных и нелинейных уравнений"

Цель занятия

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

Сценарий лабораторной работы

Указания

Метод Холецкого позволяет найти такую матрицу $$L$$, которая имеет вид$$A=LL^T,$$ где $$L$$ - нижняя треугольная матрица$$L=\left(% \begin{array}{cccc} l_{11} 0 \dots 0 \\ l_{21} l_{22} \dots 0 \\ \dots \dots \dots \dots \\ l_{n1} l_{n2} \dots l_{nn} \\ \end{array}% \right)$$ Если мы построим разложение $$A=LL^T$$, то уравнение$$Ax=f$$ может быть заменено двумя уравнениями, которые можно последовательно решить$$Ly=b,\quad L^Tx=y.$$