Материал: Спецглавы высшей математики. численные методы. Пантелеев И.Н

min ( f(х)) = f1(х)2 + f2 (x)2 +... f (х)2 + L ,

где в общем случае f(x) — вектор-функция, х — векторстолбец искомых переменных, L — некоторая константа. Синтаксис функции lsqnonlin ( ) имеет следующий вид:

х = lsqnonlin(fun,x0,eb,ub,options,PI,P2,...) [x,resnorm,residual,exitflag,output,lambda, jacobian]=lsqnonlin(...)

Он аналогичен синтаксису функции fsolve ( ), подробно описанной в разделе IV. Поэтому далее мы ограничимся только примером, демонстрирующим использование данной

функции для нахождения параметров функции

F (х,а,Ь,с) = ехр(а + bх + сх2).

Пример 7.3. Найти, используя метод наименьших квадратов, для данных, представленных в табл. VII.3, коэффициенты

аппроксимирующей функции

F (х,а,b,с) = ехр(а + bх + сх2).

Таблица VII.3

Исходные данные примера VII.3

1. Создайте файл F77.m (листинг VII.1), содержащий описание функции, возвращающей значения вектор-функции fх).

Листинг VII.1. Файл F77.m function z=F77(Coeff,vx,vy) k=l:length(vx);

z=vyexp(Coeff(1)+Coeff(2)*vx+Coeff(3)*vx.^2);

168

12

Рис. VII.4. Исходные данные и аппроксимирующая функция вида F(x,a,b,c) = exp(a + bx + cх2)

2. Выполните следующую последовательность команд:

% задание исходных данных

» vx=[0.3;0.4;l;1.4;2;4] VX =

0.3000

0.4000

1.0000

1.4000 169

2.0000

4.0000

» vy=[9.4;11.2;5;3;6;0.2] vy = 9.4000

11.2000

5.0000

3.0000

6.0000

0.2000

» z=[l 0 -1] % начальное приближение z =

1 0 -1

%вычисление коэффициентов аппроксимирующей

%функции

» Coeff = lsqnonlin('F77,z',[],[],[],vx,vy) Optimization terminated successfully: Relative function value changing by less than OPTIONS.TolFun

Coeff =

2.5696 -0.8037 0.0462

»F=inline('ехр(а+Ь*х+с*х.^2)','x','a','b','c'); % задание % аппроксимирующей функции

»X=vx(l):0.01:vx(length(vx));

%координаты абсцисс, в которых

%будут вычисляться значения

%аппроксимирующей функции

»Y=feval(F,X,Coeff(1),Coeff(2),Coeff(3)); % вычисление значений % аппроксимирующей функции

»i=l:length(vx);

»j=1:length(X);

»plot(vx(i),vy(i),'o',X(j) ,Y(j))

%визуализация исходных данных и

%аппроксимирующей функции (рис. VII.4)

170

VIII. ЧИСЛЕННЫЕ МЕТОДЫ РЕШЕНИЯ ОБЫКНОВЕННЫХ ДИФФЕРЕНЦИАЛЬНЫХ УРАВНЕНИЙ

Рассмотрим подходы и методы решения обыкновенных дифференциальных уравнений (методы Пикара, Эйлера, ЭйлераКоши, Рунге-Кутта), а также обсудим использование соответствующих функций пакета MATLAB.

8.1. Общие сведения и определения

Определение VIII.1. Дифференциальным уравнением 1-го порядка называют соотношение вида

Определение VIII.2. Дифференциальное уравнение вида

где f(x,y) — заданная функция двух переменных, называется дифференциальным уравнением первого порядка, разрешенным относительно производной.

Определение VIII.3. Решением дифференциального уравнения на интервале I называется непрерывно дифференцируемая функция у = φ(х), превращающая уравнение в тождество на I.

Для дифференциального уравнения первого порядка (VIII.1), (VIII.2), по определению получаем:

Определение VIII.4. Соотношение (VIII.3) называется решением уравнения (VIII.4) в неявной форме (или интегралом уравнения (VIII.2)), если оно определяет у как функцию от x: у = φ(x), которая есть решение уравнения

(VIII .2).

Определение VIII.5. График решения у = φ(x) уравнения (VIII.2) называется интегральной кривой данного уравнения. Определение VIII.6. Проекция графика решения на ось ординат называется фазовой кривой или траекторией дифференциального уравнения.

Определение VIII.7. Задача о нахождении решения у = φ(x) уравнения (VIII.2), удовлетворяющего начальному условию φ (хо)=уо, называется задачей Коши.

Определение VIII.8. Через каждую точку (x,y) из области определения уравнения (VIII.2) проведем прямую, тангенс угла которой к оси абсцисс равен f(x,y). Данное семейство прямых называется полем направлений, соответствующим уравнению (VIII.2) (или полем направлений функции f(x,y)).

Интегральная кривая в каждой своей точке касается поля направлений функции f(x,у).

Существование и единственность задачи Коши дифференциальных уравнений (VIII.1), (VIII.2) обеспечивается теоремой Пикара.

Теорема Пикара. Если функция f определена и непрерывна в некоторой области G, определяемой неравенствами

и удовлетворяет в этой области условию Липшица по у:

то на некотором отрезке │x - xo│≤ h, где h — положительное число, существует и притом только одно решение у = у(х) уравнения (VIII.2), удовлетворяющее начальному условию у0 = у(х0). Здесь М — константа Липшица, зависящая в общем случае от а и b. Если f(x,y) имеет в G ограниченную

172

Определение VIII.9. Дифференциальным уравнением n-го порядка называют соотношение вида

где у0 , у’0 , у0(n-1), - заданные числа, называется задачей Коши для системы дифференциальных уравнений.

Разрешив дифференциальное уравнение (VIII.8) относительно производной у(n) и выполнив следующую замену переменных:

y z1 ,

y z2 , (VIII.10)

y n 1 zn 1 ,

дифференциальное уравнение n-го порядка сводится к системе дифференциальных уравнений первого порядка:

y z1, y z2 ,

Например, уравнение второго порядка

173

□Аналитические методы решения.

□Графические методы.

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

8.2. Метод Пикара

Метод Пикара позволяет получить приближенное решение дифференциального уравнения (VIII.2) в виде функции, заданной аналитически.

Пусть в условиях теоремы существования требуется найти решения (VIII.2) с начальным условием уо = у(хо). Запишем уравнение (VIII.1) в следующем эквивалентном виде:

y0 x0

Вычислив интеграл в правой части, получим:

Очевидно, что решение интегрального уравнения (VIII.16) будет удовлетворять дифференциальному уравнению (VIII.2) и начальному условию уо = у(хо).

Действительно, при х = х0 получим:

y x0 y0 x f (x, y)dx y0 .

x0

Интегральное уравнение (VIII.16) позволяет использовать

174

метод последовательных приближений. Положим у = у0 и получим из (VIII.16) первое приближение:

y1 x y0 x f (x, y0 )dx. (VIII.17)

x0

Интеграл, стоящий в правой части (VIII.17), содержит только переменную х, после нахождения этого интеграла будет получено аналитическое выражение приближения у1 как функции переменной х. Заменим теперь в уравнении (VIII.16) у найденным значением, y1(x) и получим второе приближение:

Последовательное применение формулы (VIII.19) дает последовательность функций:

y1 x , y2 x ,..., yn x ,... (VIII.20)

Так как функция f непрерывна в области G, то она ограничена в некоторой области G' G, содержащей точку (хо,y0), т. е.

Применяя к уравнению (VIII.19) принцип сжимающих отображений, можно показать, что последовательность (VIII.20) сходится по метрике

1 , 2 max 1 x 2 x

впространстве непрерывных функций φ, определенных на сегменте x x0 d , таких, что y0 Nd. Пределx

последовательности является решением интегрального

175

уравнения (VIII.16), а, следовательно, и дифференциального уравнения (VIII.2) с начальными условиями уо = у(хо). Это означает, что к-й член последовательности (VIII.20) является приближением к точному решению уравнения (VIII.2) с определенной степенью точности.

Оценка погрешности к-ro приближения дается формулой:

где М — константа Липшица (IX.7), N — верхняя грань модуля функции f из неравенства (IX.21), а величина d для определения окрестности

│x-xo│<d

вычисляется по формуле:

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

В основе метода Эйлера лежит идея графического построения решения дифференциальных уравнений (рис. VIII.1).

Рис. VIII.1. Графическая интерпретация метода Эйлера

176

Пусть дано уравнение (VIII.2) с начальным условием у0 = у(х0).

Вместо искомой интегральной кривой на отрезке [хо ,х1]

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

(VIII.2), состоит в последовательном применении пары формул:

Метод Эйлера, как видно из рис VIII.1, имеет погрешность. Найдем локальную погрешность, присутствующую на каждом шаге, которая определяется разностью между точным значением функции и соответствующим значением касательной. Для первого шага:

y x1 y0 hf x0 , y0 y x0 h y0 hf x0 , y0