Метод Рунге-Кутта для задачи Коши 1-ого порядка

Решение задачи Коши

Метод Эйлера

Описание метода

Пусть дана задача Коши для уравнения первого порядка

y’= f(x, y)

на [a, b]

y(a)=y0

Решение разыскивается на интервале [a,b]. На этом интервале введем узлы с шагом h:

y’ = lim(), при и y’

откуда мы получаем:

yi’ = , где h=

откуда мы выражаем y(i+1):

y(i+1) = yi + h*f (xi,yi) - формула метода Эйлера, основанном на аппроксимации интегральной кривой кусочно-линейной функцией.

Задание:

a=0;

b=pi/4;

f=inline('(cos(x)).^2-y.*tan(x)');

y0=0;

n=10;

P=inline('sin(x).*cos(x)');

Программа:

function [M,E]=eiler1(a,b,f,y0,n,P)

x=linspace(a,b,n+1);

h=abs(b-a)/n;

y(1)=y0;

for i=1:n;

y(i+1)=y(i)+h*f(x(i),y(i));

end

M=[x;y]';

E=max(abs(P(x)-y));

plot(x,P(x),'',x,y,'rx')

hold on

Эйлер 2-ого порядка (модифицированный)

Теория:

Аналогичен обычному методу Эйлера, отличие состоит в том что мы находим сначала промежуточные значения x (i+1/2) и y (i+1/2).

,где L1 – касательная в начальной точке. В промежуточной точке Р(с шагом h/2) построим прямую L*, где tg = f(xi+ h/2,yi+ h/2 ∙ f(xi, yi)).

И через начальную точку построим прямую L0 параллельную L*. Из графика мы видим, что погрешность уменьшилась→то есть данный метод точнее предыдущего.

Программа:

function [M,E]=eiler2(a,b,f,y0,n,P)

x=linspace(a,b,n+1);

h=abs(b-a)/n;

y(1)=y0;

for i=1:n;

xx=x(i)+h/2;

yy=y(i)+(h/2)*f(x(i),y(i));

y(i+1)=y(i)+h*f(xx,yy);

end

M=[x;y]';

E=max(abs(P(x)-y));

plot(x,P(x),'',x,y,'rx')

hold on

Эйлер 3

Первоначально грубо находим по методу Эйлера:

y(i+1)=y(i)+h*f(x(i),y(i))

дальше используем формулу со среднеарифметическим значением производной:

y(i+1)=y(i)+(h/2)*(f(x(i),y(i))+f(x(i+1),y(i+1)))

Программа:

function [M,E]=eiler3(a,b,f,y0,n,P)

x=linspace(a,b,n+1);

h=abs(b-a)/n;

y(1)=y0;

for i=1:n;

y(i+1)=y(i)+h*f(x(i),y(i));

y(i+1)=y(i)+(h/2)*(f(x(i),y(i))+f(x(i+1),y(i+1)));

end

M=[x;y]';

E=max(abs(P(x)-y));

plot(x,P(x),'',x,y,'rx')

hold on

Оцениваем погрешность по:

[M,E1]=eiler1(a,b,f,y0,n,P);

[M,EE1]=eiler1(a,b,f,y0,2*n,P);

r1=E1/EE1

[M,E2]=eiler2(a,b,f,y0,n,P);

[M,EE2]=eiler2(a,b,f,y0,2*n,P);

r2=E2/EE2

[M,E3]=eiler3(a,b,f,y0,n,P);

[M,EE3]=eiler3(a,b,f,y0,2*n,P);

r3=E3/EE3

Метод Рунге-Кутта для задачи Коши 1-ого порядка

Рассмотрим задачу Коши

y’ = f(x,y)

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

Вычисление нового значения проходит в четыре стадии:

где h — величина шага сетки по x.

Программа:

function [R,E]=rungkut(a,b,f,y0,P,n)

x=linspace(a,b,n+1);

h=abs(b-a)/n;

y(1)=y0;

for i=1:n;

k1=h*f(x(i),y(i));

k2=h*f(x(i)+h/2,y(i)+k1/2);

k3=h*f(x(i)+h/2,y(i)+k2/2);

k4=h*f(x(i)+h,y(i)+k3);

y(i+1)=y(i)+(k1+2*k2+2*k3+k4)/6;

end

R=[x;y]';

E=max(abs(P(x)-y));

plot(x,P(x),'',x,y,'rx')

hold on

Дата добавления: 2015-12-08