В данной курсовой работе будет рассмотрен метод продолжения решения по параметру, с помощью которого можно эффективно находить корни нелинейных САУ. В работе исследуется влияние вектора начальных приближений x⁰ и заданной точности решения ε_gonна число итераций, время счета и сходимость метода. Так же дается описание программного обеспечения и тексты программ, использованные в данной работе для построения графиков сходимости метода для различных начальных значений вектора x⁰, графики ошибки.

1. Постановка задачи (математическое описание метода)

Метод продолжения решения по параметру является наиболее универсальным при решении нелинейных САУ. Пусть t - параметр, меняющийся от 0 до1. Введем в рассмотрение некоторую САУ

H (x, t) =0,

такую, что:

1) При t=0 система H (x, 0) =0 имеет решение x⁰;

2) При t=1 система H (x, 1) =0 имеет решение x^*;

3) Вектор-функция H (x, t) непрерывна по t. Тогда меняя t от 0 до 1 и решая для каждого t_iсистему H (x, t_i) =0, например, методом Ньютона, можно найти последовательно x⁰, x¹, x², …, x^*.

Так как x⁰ при t=0 известно, то всегда можно найти t₁, достаточно близкое к t₀, при котором будут выполняться условия сходимости, например, метода Ньютона. Аналогично можно обеспечить условия сходимости метода Ньютона и для t₂, t₃,…, t=1.

Вектор-функция H (x, t) может быть выбрана различными способами. Рассмотрим три распространенных варианта:

1) H (x, t) =F (x) + (t-1) *F (x⁰) =0

При t=0 получаем: F (x⁰) - F (x⁰) =0, т.е. условие 1) выполнено.

При t=1 F (x^*) - (1-1) * F (x⁰) =F (x^*) =0. И, наконец, вектор-функция H (x, t) непрерывна по t.

2) H (x, t) =t*F (x).

Условия 1) - 3) соблюдаются и для этой вектор-функции.

Идея метода состоит в следующем. Полагаем t₁₌∆t и решаем систему H (x, t₁) =0 при выбранном x⁰. Получаем x^t₁. Далее, берем его в качестве начального приближения и решаем при новом t₂=t₁+∆t систему H (x, t₂) =0, получаем x^t₂ и так далее до тех пор, пока не будет достигнута заданная точность. Нелинейные системы H (x, t_i) =0 на каждом шаге по t решаются, например, методом Ньютона, который обычно сходится, так как x^t_i-1 и x^t_i лежат близко друг к другу. Если несмотря на это решение x^t_i не получается за 6-7 итераций, ∆t уменьшается и система H (x, t_i) =0 решается снова.

Последовательность шагов реализации алгоритма состоит в следующем:

Шаг 1. Формирование системы H (x, t) =0.

Шаг 2. Выбор начального приближения x⁰, (например, x⁰=0) и точности решения ε_gon.

Шаг 3. Полагаем i=1.

Шаг 4. Вычисляем t_i=t_i-1+∆t (обычно вначале берут ∆t=0,1)

Шаг 5. Решаем систему H (x, t_i) =0. Получаем вектор x^t_i. При этом считаем число итераций m. Если m>10, значит метод Ньютона уже не сойдется, так как x^t_i-1 и x^t_i слишком далеки друг от друга. Тогда надо уменьшить ∆t в два раза и вернуться к шагу 4. Будем считать, что x^t_i найдено.

Шаг 6. Проверяем, достигли ли мы заданной точности. Например, используя первый способ,

|| x^t_i-x^t_i-1|| ≤ ε_gon.

Если последнее условие не соблюдается, то переходим к шагу 4. Иначе считаем, что x^*=x^t_i и расчеты закончены.

2. Описание программного обеспечения

2.1 Общие сведения и требования к ПО и описание логической структуры

ПО состоит из следующих файлов: mpr. m, prog. m, funf. m, funj. m. Программы, реализующие метод, разработаны в среде МаtLab, предназначенной для выполнения математических операций. Программа состоит из программы-функции mpr. m, которая описывает метод, программы с данными - основная программа prog. mи двух подпрограмм-функций funf. m - для нахождения корней системы уравнений; funj. m - для нахождения матрицы Якоби. Рассмотрим их подробнее.

Функциональное назначение

Программа предназначена для решения систем нелинейных алгебраических уравнений в среде МаtLab методом продолжения решения по параметру.

Используемые переменные:

t - время выполнения итерационного процесса;

x - вектор начального приближения к решению;

n - размерность вектора;

m - номер итерационного процесса;

it - счетчик итераций.

Входные параметры:

funf - формальное имя программы, которое дает возможность вычислить корни нелинейных САУ.

funj - формальное имя программы, которое дает возможность вычислить матрицу Якоби.

x0 - начальное приближение собственного вектора;

dt - приращение времени;

edop - заданная допустимая ошибка;

trace - установка режима вывода на экран;

Выходные параметры:

tout - выходное значение времени;

xout - конечное значение x;

dxout - конечное значение вектора ошибки.

Тексты программ:

Mpr. m

function [xout,dxout,tout] =mpr (funf,funj,x0,dt,edop,trace)

t=dt; x=x0; tout=t; xout=x0'; n=size (x0);

dxout=zeros (1,n); m=0; it=0;

f0=feval (funf,x0);

while (t<=1)

ndx=1;

nh=1;

nv= [ndx; nh];

while (max (nv) >edop)

J=feval (funj,x0);

F=feval (funf,x0);

h= (-F) *t;

dx=J\h;

x=x+dx;

m=m+1;

ndx=norm (dx);

nh=norm (h);

nv= [ndx; nh];

if (m > 10)

t=t-dt;

dt=dt/2;

t=t+dt;

x=x0;

m=0;

end;

x0=x;

tout= [tout; t];

xout= [xout; x'];

dxout= [dxout; dx'];

if (m < 4)

dt=dt*2;

end;

t=dt+t;

it=it+1;

end;

disp ('it='); %количество итераций

disp (it);

disp ('t='); %время выполнения итерационного процесса

disp (t);

pause;

xout - конечное значение x;

dxout - конечное значение вектора ошибки.

m - номер итерации.

Prog. m

trace=1;

dt=0.1;

x0=0;

edop=0.1;

[xout,dxout,m] = mpr ('funf','funj',x0,dt,edop,trace);

plot (m,xout); %График значений x

pause;

plot (m,dxout); %График ошибки

pause;

Funf. m

function [f] =funf (x)

f= [0.0001*exp (30*x) +x-6];

end

Funj. m

function [j] =funj (x)

j= [30*0.0001*exp (30*x) +1];

end

3. Описание тестовых задач

Исследуем влияние вектора начальных приближений на время счета, число итераций и сходимость метода.

Начальное приближение x⁰=0, заданная точность edop=0.1, dt=0.1.