Численные методы решения ОДУ Опции

[Veina]

Я написала программу на Паскале для формулы y_n+1=y_n+h*f_n
Вот код:

program eiler;

var x,y,a,b,h:real;

function f(x,y:real):real;

begin f:= y*y*(8-20*x); end;

begin

writeln('vvedite y, a, b, h');

readln(y,a,b,h); x:=a;

repeat

writeln(x:0:1,' ',y:0:3);

y:=y+h*f(x,y);

x:=x+h;

until not (x<b);

readln;

end.

А как написать переделать для формулы
y_n+1=y_n+h*f_n+1
и для формулы y_n+1=y_n+(h/2)*(3f_n-f_n-1)
Подскажите, пожалуйста, кто-нибудь. Заранее спасибо

[Lapp]

Для f(x[n], y[n]) - да, было.
А для f(x[n+1], y[n]) - не уверен. Я такого не помню.

А зачем помнить? Можно порассуждать )).
1. Сначала представим себе, что f не зависит от y. Если мы аппроксимируем y прямой y'*x на участке от x_n до x_n+1, то неважно, где брать значение y' - слева или справа. То есть обе формулы имеют право на жизнь и примерно одинаковую точность. Было бы лучше, если бы мы брали значение f в центре. Логично также взять среднее ар. значений на концах (или другое взвешенное значение). Есть также методы, которые предполагают наличие двух различных сеток для x и y. Добавление зависимости от y дела не меняет (в определенных предположениях о функции, ессно). Так или иначе, сходимость должна быть.

2. Я вообще не вижу смысла вникать в выч.мат. в данном случае. Человек принес с собой формулы и просит помочь с программированием их. Если ты считаешь, что есть формулы лучше - ну, упомяни об этом... Но это, мне кажется, не секрет )).

3. С тем, что для реализации последней формулы нужно как-то (по-другому) рассчитывать первое значение, ты, пожалуй, прав.. Я как-то не сразу врубился, о чем речь, извини. Об этом я сейчас еще напишу (ниже).

Тогда для формулы y_n+1=y_n+(h/2)*(3f_n-f_n-1)
будет вот так:
<...>
Правильно? проверьте, пожалуйста

Нет, не правильно.
Во-первых, ты воплотила формулу:
y_n+1=y_n+(h/2)*(3f_n+1-f_n)
- а это СОВСЕМ не то..
Марин, тут есть две возможности.
1. Если нам разрешено использовать функцию f за пределами отрезка [a,b], то все просто:

program eiler_2;

var
  x,y,a,b,h: real;

function f(x,y: real): real;
begin
  f:= y*y*(8-20*x);
end;

begin
  writeln('vvedite y, a, b, h');
  readln(y,a,b,h);
  x:= a;
  repeat
    writeln(x:0:1,' ',y:0:3);
    y:= y+(h/2)*(3*f(x,y)-f(x-h,y));    // исправлена формула
    x:= x+h;
  until not (x<b);
  readln;
end.

- и все.

2. Если за пределы отрезка [a,b] по x вылезать нельзя, то для реализации этой формулы нужно сначала подсчитать y₁ (отсчет от 0) каким-то другим методом, как уже говорил TarasBer. Каким именно - не сказано, так что ты можешь выбрать любой (первый, например). Точность решения от этого сильно зависеть не будет все равно.
Примерно вот так:

program eiler_3;

var
  x,y,a,b,h: real;

function f(x,y: real): real;
begin
  f:= y*y*(8-20*x);
end;

begin
  writeln('vvedite y, a, b, h');
  readln(y,a,b,h);

  writeln(a:0:1,' ',y:0:3);   // печать перед первой итерацией
  y:= y+h*f(x,y);              // первая итерация (первый метод)

  x:= a+h;                    // цикл начинаем со второй итерации
  repeat
    writeln(x:0:1,' ',y:0:3);
    y:= y+h/2*(3*f(x,y)-f(x-h,y));    // честно пишем формулу
    x:= x+h;
  until x>=b;
  readln;
end.

Понятно? Спрашивай, что неясно.