Технический форум - Подынтегральная функция методом Симпсона

Технический форум (http://www.tehnari.ru/)

- Delphi, Kylix and Pascal (http://www.tehnari.ru/f43/)

Подынтегральная функция методом Симпсона

Помогите с задачей и объясните ее доступным для новичка языком:tehnari_ru_674:

Задача: вычисление определённых интегралов и табулирование первообразных функций

Все что мне известно:

Подынтегральная функция - 1/x*ln*x
Промежуток интеграла - [2;3]
Кол-во частей разбивки - 36
Вычисление значения первообразной h - 0,2
Точность вычисления значения первообразной - 0,001
Точное значение первообразной - 2,3026(ln ln x - ln ln 2)

Да легко, сейчас нарисую.

Значит, так.
Прежде всего — методичка составлена бестолково и ошибками. В частности, откуда взялся множитель 2,3026 — совершенно не представляю. А так — пожалуйста:

Код:

Const

 a=2.0;

 Eps=0.001;

 h=0.2;



Var

 x:Real;

 j:Byte;



Function F(p:Real):Real;

begin

 F:=1.0/(p*ln(p));

end;



Function Antiderivative(p:Real):Real;

begin

 Antiderivative:=Ln(Ln(p));

end;





Function Simpson(p:Real):Real;

var

 i,N:Integer;

 Sum1,Sum2,t:Real;

begin

 N:=2;

 t:=(p-a)/(2*N);

 Sum1:=(F(a)+F(p))/3*t;

 for i:=1 to N do

  Sum1:=Sum1+4.0*t/3.0*F(a+t*(2*i-1));

 for i:=2 to N do

  Sum1:=Sum1+2.0*t/3.0*F(a+t*(2*i-2));

 Sum2:=Sum1;

 Repeat

  Sum1:=Sum2;

  N:=N*2;

  t:=(p-a)/(2*N);

  Sum2:=(F(a)+F(p))/3*t;

  for i:=1 to N do

   Sum2:=Sum2+4.0*t/3.0*F(a+t*(2*i-1));

  for i:=2 to N do

   Sum2:=Sum2+2.0*t/3.0*F(a+t*(2*i-2));

 Until ABS(Sum2-Sum1)<Eps;

 Simpson:=Sum2;

end;



Begin

 WriteLn('   x     Simpson    Accurate value');

 for j:=1 to 5 do

  begin

   x:=a+h*j;

   WriteLn(x:5:1, '   ',Simpson(x):0:8, '   ',(Antiderivative(x)-Antiderivative(a)):0:8);

  end;

 ReadLn

End.

Спасибо вам)) Я как дома буду проверю
А как вам мой вариант решения задачи, помог препод:

Цитата:

Сообщение от Alessandro4 (Сообщение 2632167)

А как вам мой вариант решения задачи

Извините, но вижу лишь обрывок весьма, мягко говоря, странного кода. Симпсоном тут и не пахнет.

Цитата:

Сообщение от Alessandro4 (Сообщение 2632167)

Я как дома буду проверю

А что проверять? Я ж выложил результат работы программы.

Большое вам спасибо!
Я больше имел в виду что изучу ее, чтоб суметь грамотно объяснить преподавателю

Цитата:

Сообщение от Alessandro4 (Сообщение 2632194)

Я больше имел в виду что изучу ее, чтоб суметь грамотно объяснить преподавателю

Удачи! И обращайтесь, если потребуются с моей стороны какие-то объяснения.

Ладно, пытался я значит сам разобраться, но что то не очень, выходит:tehnari_ru_325:
Объясните, пожалуйста, ваши блоки и сам метод)

Цитата:

Сообщение от Alessandro4 (Сообщение 2632299)

Объясните, пожалуйста, ваши блоки и сам метод)

Пожалуйста.
Интеграл берется аналитически, и первообразная (кстати, в методичке, как обычно бывает, это слово употреблено безграмотно, на самом деле первообразная, она же "неопределенный интеграл", это функция, дифференцирование которой даёт подынтегральное выражение, поэтому первообразная определяется с точностью до произвольной константы) есть
Ln(Ln(x))+С.
Отсюда определенный интеграл с переменным верхним пределом (который эти умники называют "первообразной") есть
Ln(Ln(x)) - Ln(Ln(a)).
Уж откуда эти ... взяли числовой коэффициент — не ведаю. Впрочем, догадываюсь: скорее всего, подынтегральная функция предполагалась не 1/(xLn(x)), как у них в таблице, а 1/(xLg(x)), т.е. с десятичным, а не натуральным логарифмом. Тогда действительно возникает коэффициент Ln(10). Идиоты!
Всё это считается в программе через функцию Antiderivative.

Теперь, собственно, о численном интегрировании.

1. Вспомним формулу Симпсона для вычисления интеграла на промежутке (a,b):
а) разбиваем промежуток на 2N частей (тут важна чётность).
б) вводим
h = (b - a)/(2N)
Xj = a + h*j, j = 0, 1, ... 2N
Yj = f(Xj), где f — подынтегральная функция.
тогда интеграл есть
(h/3)*(Y0 + 4*Y1 + 2*Y2 + 4*Y3 + ... + Y2N)

2. Тут опять "уважаемые" преподаватели бьют бедных студиозусов ключом (разводным или шведским) по голове. Видите ли, можно либо зафиксировать число разбиений интервала, тогда точность — уж как получится, либо задать точность и дробить интервал до тех пор, пока мы этой точности не достигнем. А то и другое сразу, как предписывают эти корифеи, НЕВОЗМОЖНО. Поэтому я выбрал второй вариант — вычисление интеграла с заданной точностью. Вычисление делается в функции Simpson.

3. В этой функции организован следующий цикл:
а) Вычисляется значение интеграла при 2N = 2. Обозначим это значение Sum1.
б) Удваиваем N и повторяем вычисление. Получаем величину Sum2.
в) Если |Sum1-Sum2|< ε, то на этом заканчиваем, в противном случае величине Sum1 присваиваем значение Sum2, снова вдвое дробим интервал и повторяем вычисление, в результате которого получаем новое значение Sum2, которое сравниваем с Sum1.
г) Делаем так до тех пор, пока не выполнится условие малости разницы результатов двух итераций.

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

Код:

 begin

1  N:=2;

2  t:=(p-a)/(2*N);

3  Sum1:=(F(a)+F(p))/3*t;

4  for i:=1 to N do

5   Sum1:=Sum1+4.0*t/3.0*F(a+t*(2*i-1));

6  for i:=2 to N do

7   Sum1:=Sum1+2.0*t/3.0*F(a+t*(2*i-2));

8  Sum2:=Sum1;

9  Repeat

10   Sum1:=Sum2;

11   N:=N*2;

12  t:=(p-a)/(2*N);

13   Sum2:=(F(a)+F(p))/3*t;

14  for i:=1 to N do

15   Sum2:=Sum2+4.0*t/3.0*F(a+t*(2*i-1));

16   for i:=2 to N do

17   Sum2:=Sum2+2.0*t/3.0*F(a+t*(2*i-2));

18 Until ABS(Sum2-Sum1)<Eps;

19  Simpson:=Sum2;

20 end;

Спасибо за разъяснение, но многое мне, конечно, не понятно
Кстати и препод сказал, что вы ошиблись хех, да да... а именно в строчке - 7, опираясь на формулу сказал, что "-2" не должно там быть и попросил меня понять как работает весь этот цикл > объяснить ему > сдать, ибо мое объяснение его не особо удовлетворило. Объясните, пожалуйста, все строчки, если вас не затруднит, ибо в методе симпсона и в самих интегралах я не очень, желательно разжевав) Еще раз огромное спасибо)