Интегральные уравнения

@SW Developer · Регистрация: 10.01.2016

Студворк — интернет-сервис помощи студентам

Здравствуйте!

Нужно решить интегральные уравнения в среде Паскаль. Методами : Правых прямоугольников; Левых прямоугольников; Средних прямоугольников; Трапеций; Методом Симпсона; Методом квадратур Гауса – Лежандра; Методом Монте-Карло.

У кого есть желание поразмыслить, буду рад Вашей помощи.

Админ! Уравнение слишком большое, чтобы набирать его в редакторе формул, так что sorry.

@SW Developer · 24.03.2017, 11:48 **[ТС]**

Видно я чего-то не понимаю.
Мне просто интересно узнать, почему одним методом (правые прямоугольники) я получаю конечный результат в виде числового значения, а другим методом (левых прямоугольников) у меня нет этого значения, а получается inf. Напрашивается вывод, что один из методов не работает для данного выражения. Или я ошибаюсь?

@bormant · 24.03.2017, 11:56

Сообщение от SW Developer

Или я ошибаюсь?

Да, ошибаетесь.
Все эти методы являются методами приближенного вычисления определенного интеграла.
Все они имеют погрешность. Оценка погрешности производится на основании производных первого и следующего порядков.
При вычислении методом правых прямоугольников вы получили число с погрешностью, многократно превышающей это полученное число. Пригодность такого результата -- нулевая.
Вспомните график y=1/x. В одном случае у вас прямоугольник опирается на точку на графике, а в другом пытается опереться на inf (при 1/0).

@SW Developer · 24.03.2017, 15:25 **[ТС]**

Вот так работает:

Pascal
1
2
3
4
5
6
7
8
function f5(x: Real): Real;
         begin
              x:=a*x;
              if x<0 then x:=1/Exp(x)-1;
              if x=0 then exit
                     else x:=Exp(x)-1;
              f5:=1/x;
         end;

@bormant · 24.03.2017, 17:01

Сообщение от SW Developer

Вот так работает

... вот так работает неправильно -- для бесконечно большой величины дает 0.
Да, так оно не падает от деления на 0. Но если на входе интервал, приводящий к делению на 0, то проблема именно тут, в самом интервале, а не внутри функции.
Если уж на то пошло, и отдельные функции проверки на ОДЗ аргумента делать лень, можно сделать так:

Pascal
1
2
3
4
5
6
7
8
9
10
11
12
var fErr: Boolean;
const mxDbl = 1.7e308;
 
function f5(x: Real): Real;
begin
  x:=a*x;
  if x<0 then x:=1/Exp(x)-1 else x:=Exp(x)-1;
  if x=0 then begin
    f5:=mxDbl; fErr:=True; Exit;
  end;
  f5:=1/x;
end;

Добавлено через 3 минуты

Pascal
1
2
3
fErr:=False;
WriteLn(IntRL(...));
if fErr then WriteLn('Внимание! При вычислении значения возникло деление на 0!');

@SW Developer · 24.03.2017, 18:26 **[ТС]**

Интересно! Спасибо. Завтра проверю.

Добавлено через 7 минут
Вообще по условию, интеграл данного выражения вычисляется на отрезке [0, 7]. Если изначально известно, что функция при х=0 стремится к бесконечности, то для чего все это? Какая-то бессмыслица! Ведь и вся площадь, ограниченная графиком и осью Х стремится к бесконечности, следовательно и интеграл стремится к бесконечности.

@Cyborg Drone · 24.03.2017, 19:15

SW Developer, несколько неправильное у Вас представление о несобственных интегралах. Для многих функций вполне может быть, что в каких-то точках значение функции стремится к бесконечности, но, тем не менее, определённый интеграл вполне себе конечный. К примеру, возьмите функцию y=1/x. Если ни один из пределов интегрирования не равен 0, то определённый интеграл от 1/x имеет конечное значение.

Естественно, такие рассуждения справедливы для аналитического, а не для численного решения.

Отмечу, что интеграл

$https://www.cyberforum.ru/cgi-bin/latex.cgi? \int_{0}^{7}\frac{dx}{e^{ax}-1}\ \rightarrow \ \infty $

действительно расходится (стремится к бесконечности).

Как пример,

$https://www.cyberforum.ru/cgi-bin/latex.cgi? \int_{0}^{7}\frac{dx}{\left( e^{ax}-1\right)^x}\ \approx \ 1.31591 $

несмотря на то, что в точке x=1 подынтегральная функция стремится к ±∞.

@SW Developer · 24.03.2017, 20:56 **[ТС]**

Я так понимаю, что во втором выражении в знаменателе нуль в нулевой степени при х=0?

@Cyborg Drone · 24.03.2017, 23:48

Да. И, кроме того, при x=1 знаменатель обращается в 0, делить на 0 как-то нельзя, и функция при x=1 не определена (в данном случае ещё говорят, что значение функции стремится к бесконечности). Несмотря на всё это, определённый интеграл (правильнее сказать - определённый несобственный интеграл второго рода) сходится, иными словами, вопреки всему "очевидному", площадь между кривой и осью x конечна, и равна выше написано чему.

@SW Developer · 25.03.2017, 00:21 **[ТС]**

Т.е. такой вариант кода программы можно считать правильным

Pascal
1
2
3
4
5
6
7
8
function f5(x: Real): Real;
   begin
         x:=a*x;
         if x<0 then x:=1/Exp(x)-1;
         if x=0 then exit
                    else x:=Exp(x)-1;
         f5:=1/x;
   end;

Добавлено через 7 минут
...?

@Cyborg Drone · 25.03.2017, 00:58

Тут такое дело... Можно, конечно, считать это правильным для вычисления функции... Наверное, всё-таки вот так:

Pascal
1
2
3
4
5
6
7
8
9
function f5(a, x: Real): Real;
begin
  if (x <> 0) and (a <> 0) 
    then f5 := 1 / (exp(a * x) - 1);
    else begin
      write('Ошибка: деление на 0');
      halt
    end
end;

Короче, "a" не должно быть равно 0, и пределы интегрирования не должны включать x=0, поскольку в этих случях интеграл ну никак не сходится. То есть, при a≠0 интеграл, приведённый в задании, или, к примеру (он же, но другие пределы нитегрирования):

$https://www.cyberforum.ru/cgi-bin/latex.cgi? \int_{-1}^{1}\frac{dx}{e^{ax}-1} $

вычислять смысл умер, потому что они не сходятся, в данном случае, из-за неустранимого разрыва второго рода в точке x=0, а, к примеру,

$https://www.cyberforum.ru/cgi-bin/latex.cgi? \int_{1}^{10}\frac{dx}{e^{ax}-1} $

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

@SW Developer · 25.03.2017, 01:36 **[ТС]**

И как быть, если в задании отрезок интегрирования [0, 7], а используя метод левых прямоугольников мы первым делом определяем Y при x=0?

Добавлено через 3 минуты
С параметром а проблем не будет, если в основном блоке программы прописать условие а<> 0.

Добавлено через 31 секунду
При вводе этого параметра с клавиатуры.

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

@Cyborg Drone · 25.03.2017, 09:12

При любом методе результат не будет корректным всё равно. Забудьте об вычислении этого интеграла, если x=0 входит в интервал интегрирования, это я Вам как краевед говорю. Просто при методе правых прямоугольников на самом деле будет вычислен интеграл

$https://www.cyberforum.ru/cgi-bin/latex.cgi? \int_{0+h}^{7}\frac{dx}{e^{ax}-1} $

и величина его будет всецело зависеть от шага разбиения h. И то, что получится, никак нельзя считать результатом. Уменьшил h в несколько раз (увеличил количество шагов интегрирования) - и значение интеграла существенно выросло, но, если бы интеграл сходился, значение бы существенно не изменилось. Если попытаться найти этот интеграл с заданной точностью, используя метод уменьшения шага, то, что этот интеграл расходится, сразу же выползет наружу. Для примера, воспользуемся критерием оценки точности вычислений, использующимся обычно для методов (всех) прямоугольников и метода трапеций

$https://www.cyberforum.ru/cgi-bin/latex.cgi? \frac{\,\left|I_n-I_{2n} \right|\,}{3}<\varepsilon $

где I_n - интеграл, вычисленный при (текущем) разбиении интервала интегрирования на n частей, I_2n - интеграл, вычисленный при разбиении интервала интегрирования на 2n частей, ε - желаемая точность вычислений.

Вот сейчас выяснится, что при увеличении количества разбиений указанная дробь будет увеличиваться, а не уменьшаться, и даже точность вычислений ε<0.1 (это, скажем так, плюс-минус лапоть) не будет достигнута никогда.

Вычислим Ваш интеграл методом правых прямоугольников (для простоты примем a=1), и посмотрим точность при увеличении числа разбиений (сама величина интеграла не интересна, её даже выводить не станем):

Pascal
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
function f(x: real): real;
begin
  f := 1 / (exp(x) - 1)
end;
 
function int_r(x0, x1: real; n: integer): real;
var i: integer;
    h, s: real;
begin
  s := 0;
  h := (x1 - x0) / n;
  for i := 1 to n do s := s + f(x0 + i * h);
  int_r := s * h
end;
 
const x0 = 0;
      x1 = 7;
var I, I2: real;
    n: integer;
    a: string;
begin
  n := 2;
  I2 := int_r(x0, x1, n);
  repeat
    I := I2;
    n := 2 * n;
    I2 := int_r(x0, x1, n);
    write('n = ', n:9, ', eps = ', abs(I2 - I) / 3:0:15, ' ->');
    readln(a)
  until a = 'f';
  readln
end.

Итак, что получилось:

Нажмите, чтобы посмотреть протокол работы программы

Code
n =         4, eps = 0.107069283501296 ->
n =         8, eps = 0.163406225763953 ->
n =        16, eps = 0.195869867067778 ->
n =        32, eps = 0.213122681164054 ->
n =        64, eps = 0.222001832559446 ->
n =       128, eps = 0.226504446558399 ->
n =       256, eps = 0.228771504024013 ->
n =       512, eps = 0.229908969805996 ->
n =      1024, eps = 0.230478686923845 ->
n =      2048, eps = 0.230763791537269 ->
n =      4096, eps = 0.230906405357497 ->
n =      8192, eps = 0.230977727645918 ->
n =     16384, eps = 0.231013392634794 ->
n =     32768, eps = 0.231031226090292 ->
n =     65536, eps = 0.231040143058441 ->
n =    131072, eps = 0.231044601602443 ->
n =    262144, eps = 0.231046830889612 ->
n =    524288, eps = 0.231047945536755 ->
n =   1048576, eps = 0.231048502861629 ->
n =   2097152, eps = 0.231048781523483 ->
n =   4194304, eps = 0.231048920856109 ->
n =   8388608, eps = 0.231048990519747 ->
n =  16777216, eps = 0.231049025356252 ->
n =  33554432, eps = 0.231049042766964 ->
n =  67108864, eps = 0.231049051477989 ->
n = 134217728, eps = 0.231049055833297 ->
n = 268435456, eps = 0.231049058016822 ->f

Как видите, оценка точности вычислений и не думает стремиться к 0, а, наоборот, увеличивается.

Заменим подынтегральную функцию на что-нибудь "человеческое", например,

$https://www.cyberforum.ru/cgi-bin/latex.cgi? \int_{0}^{7}x^2dx $

Для этого в программе заменим:

Pascal
3
  f := x * x

Получим:

Нажмите, чтобы посмотреть протокол работы программы

Code
n =         4, eps = 17.864583333333333 ->
n =         8, eps = 8.039062500000000 ->
n =        16, eps = 3.796223958333333 ->
n =        32, eps = 1.842285156250000 ->
n =        64, eps = 0.907185872395833 ->
n =       128, eps = 0.450103759765625 ->
n =       256, eps = 0.224179585774740 ->
n =       512, eps = 0.111871719360352 ->
n =      1024, eps = 0.055881341298421 ->
n =      2048, eps = 0.027927041053772 ->
n =      4096, eps = 0.013960113128026 ->
n =      8192, eps = 0.006979204714298 ->
n =     16384, eps = 0.003489389394720 ->
n =     32768, eps = 0.001744641456753 ->
n =     65536, eps = 0.000872307418225 ->
n =    131072, eps = 0.000436150407921 ->
n =    262144, eps = 0.000218074349344 ->
n =    524288, eps = 0.000109036939458 ->
n =   1048576, eps = 0.000054518439242 ->
n =   2097152, eps = 0.000027259201299 ->
n =   4194304, eps = 0.000013629602895 ->
n =   8388608, eps = 0.000006814800434 ->
n =  16777216, eps = 0.000003407401692 ->
n =  33554432, eps = 0.000001703693949 ->
n =  67108864, eps = 0.000000851859132 ->
n = 134217728, eps = 0.000000425934426 ->
n = 268435456, eps = 0.000000212972774 ->f

А здесь - всё нормально. Несмотря на то, что оценка точности вычислений сначала ни в какую эпсилон не лезет, она при увеличении n устремляется к 0. И заданный плюс-минус лапоть ε<0.1 был достигнут при десятом удвоении количества разбиений n.

@SW Developer · 25.03.2017, 13:24 **[ТС]**

bormant, объясните, пожалуйста:

Сообщение от bormant

const mxDbl = 1.7e308;

Добавлено через 2 минуты
Я так понимаю - интеграл стремится в бесконечность?

@bormant · 25.03.2017, 13:29

При делении 1 на бесконечно малое -- так и есть.

@SW Developer · 26.03.2017, 16:08 **[ТС]**

А алгоритм по методу квадратуры Гаусса–Лежандра есть у кого-нибудь?

@SW Developer · 29.03.2017, 12:41 **[ТС]**

Вроде как написал, но что-то сомнения закрались, уж ОООчень большая разница в вычислениях:

Pascal
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
//Метод квадратур Гаусса-Лежандра (2-х точечный)
function IntGauss(a, b: Real; dx: Word; f: TFunc): Real;
         var x0, x1, x2, x11, x12, x13, x14: Real; i: byte;
         begin
         for i:=1 to dx do begin
              x0:=0;
              x11:= (b+a)/2;
              x12:= (b-a)/2*sqrt(3);
              x13:= x11 - x12;
              x14:= x11 + x12;
              x1:= f(x13);
              x2:= f(x14);
              x0:= x0 + x1+x2;
              end;
              IntGauss:= (b-a)/2 * x0;
         end;

так сказать черновой вариант.

@SW Developer · 29.03.2017, 12:45 **[ТС]**

Формулу использовал вот эту:

$https://www.cyberforum.ru/cgi-bin/latex.cgi? I\approx \frac{b-a}{2}\left(f\left(\frac{a+b}{2}-\frac{b-a}{2\sqrt{3}} \right)+f\left(\frac{a+b}{2}+\frac{b-a}{2\sqrt{3}} \right)\right) $

@SW Developer · 29.03.2017, 13:42 **[ТС]**

Цикл не нужен!

Pascal
1
2
3
4
5
6
7
8
9
10
11
12
function IntGauss(a, b: Real; f: TFunc): Real;
         var x1, x2, x11, x12, x13, x14: Real; i: byte;
         begin
   
              x11:= (b+a)/2;
              x12:= (b-a)/2*sqrt(3);
              x13:= x11 - x12;
              x14:= x11 + x12;
              x1:= f(x13);
              x2:= f(x14);
              IntGauss:= (b-a)/2 * (x1+x2);
         end;

Добавлено через 3 минуты
Но, результат все равно в 4 раза больше, чем при вычислении другими методами.

Добавлено через 39 минут
Получилось!

Pascal
1
2
3
4
5
6
7
8
9
10
11
12
13
function IntGauss(a, b: Real; f: TFunc): Real;
         const x0 = 0.577350269;
         var   h, x1, x2, x11, x12, x13, x14: Real;
         begin
              h:=(b-a)/2;
              x11:= (b+a)/2;
              x12:= (b-a)/2*x0;
              x13:= x11 - x12;
              x14:= x11 + x12;
              x1:= f(x13);
              x2:= f(x14);
              IntGauss:= h * (x1+x2);
         end;

Новые блоги и статьи Все статьи Все блоги /
интеграция AnyLogic с самописным REST API и переход на Odoo anaschu 03.07.2026 Успешная интеграция AnyLogic с самописным REST API и переход на промышленную Odoo WMS Сегодня проделал огромный путь от простой симуляции физических процессов до построения полноценной. . .	Поиск всех путей на ориентированном графе. Linux dcc0 02.07.2026 Переработка старого кода из моей статьи. Через несколько переработок от PHP кода к C89 (надеюсь, 89). Но довольно запутанно получилось. Код для Linux. Но если убрать time и то, что с ним. . .	Сам себя обучал rest api anaschu 02.07.2026 Педагогический лайфхак: Почему чистый REST API для ученика намного круче, чем готовые библиотеки Когда мы отказались от капризного JAR-файла AnyLogic и переписали код на стандартный HttpClient,. . .	rest api anylogic - выполнение модели на своём русском сайте anaschu 02.07.2026 Как подружиться с AnyLogic Cloud API, победить провайдеров и развернуться Java-бэкенд в Docker на бесплатном хостинге: Двухдневный лог борьбы Всем привет! Хочу поделиться свежим (и довольно. . .
Где деньги лежат kumehtar 02.07.2026 Это - японская подводная лодка I-52 (тип C2, кодовое имя Momi) вышла из Японии в марте 1944 года с миссией в оккупированную немцами Францию (Лорьян). Это была одна из «Янаги»-миссий по обмену. . .	Krabik для WoW 3.3.5a, многоязычный AmbA 02.07.2026 Допилил бота, думаю что окончательно. Изменения: - добавлена многоязычность - добавлено снятие скриншотов - добавлено поддержание бафов хождения по воде (для жреца, дк и шамана) - и так, по. . .	Алиса нашла кучу ошибок компиляции и запуска в проекте, который без проблем компилировался и запускался))) anaschu 30.06.2026 Я пока посмеюся, но завтра проверю. А вообще интерсно. Дал алисе файл, в котором точно нет ошибок компиляции и запуска, и попросил их найти. Нашла кучу))) Критические ошибки, мешающие компиляции и. . .	сукцессия 16. Общий обзор, в основном что бы другие ии поняли anaschu 29.06.2026 # Передаточный документ: модель микоризной сукцессии (для нового чата) Этот документ предназначен для того, чтобы новый чат Claude мог продолжить работу без необходимости заново разбираться в. . .

@SW Developer 97 / 93 / 81 Регистрация: 10.01.2016 Сообщений: 663 Записей в блоге: 13

	Интегральные уравнения 17.03.2017, 14:24. Показов 4135. Ответов 37 Метки Интегралы (Все метки) Здравствуйте! Нужно решить интегральные уравнения в среде Паскаль. Методами : Правых прямоугольников; Левых прямоугольников; Средних прямоугольников; Трапеций; Методом Симпсона; Методом квадратур Гауса – Лежандра; Методом Монте-Карло. У кого есть желание поразмыслить, буду рад Вашей помощи. Админ! Уравнение слишком большое, чтобы набирать его в редакторе формул, так что sorry. Миниатюры 0

@SW Developer 97 / 93 / 81 Регистрация: 10.01.2016 Сообщений: 663 Записей в блоге: 13
	24.03.2017, 11:48 [ТС]
	Видно я чего-то не понимаю. Мне просто интересно узнать, почему одним методом (правые прямоугольники) я получаю конечный результат в виде числового значения, а другим методом (левых прямоугольников) у меня нет этого значения, а получается inf. Напрашивается вывод, что один из методов не работает для данного выражения. Или я ошибаюсь? 0

@SW Developer 97 / 93 / 81 Регистрация: 10.01.2016 Сообщений: 663 Записей в блоге: 13
	24.03.2017, 18:26 [ТС]
	Интересно! Спасибо. Завтра проверю. Добавлено через 7 минут Вообще по условию, интеграл данного выражения вычисляется на отрезке [0, 7]. Если изначально известно, что функция при х=0 стремится к бесконечности, то для чего все это? Какая-то бессмыслица! Ведь и вся площадь, ограниченная графиком и осью Х стремится к бесконечности, следовательно и интеграл стремится к бесконечности. 0

@Cyborg Drone Модератор 10465 / 5761 / 3410 Регистрация: 17.08.2012 Сообщений: 17,515
	24.03.2017, 19:15
	SW Developer, несколько неправильное у Вас представление о несобственных интегралах. Для многих функций вполне может быть, что в каких-то точках значение функции стремится к бесконечности, но, тем не менее, определённый интеграл вполне себе конечный. К примеру, возьмите функцию y=1/x. Если ни один из пределов интегрирования не равен 0, то определённый интеграл от 1/x имеет конечное значение. Естественно, такие рассуждения справедливы для аналитического, а не для численного решения. Отмечу, что интеграл $https://www.cyberforum.ru/cgi-bin/latex.cgi?<br /> \int_{0}^{7}\frac{dx}{e^{ax}-1}\ \rightarrow \ \infty<br />$ действительно расходится (стремится к бесконечности). Как пример, $https://www.cyberforum.ru/cgi-bin/latex.cgi?<br /> \int_{0}^{7}\frac{dx}{\left( e^{ax}-1\right)^x}\ \approx \ 1.31591<br />$ несмотря на то, что в точке x=1 подынтегральная функция стремится к ±∞. 0

@SW Developer 97 / 93 / 81 Регистрация: 10.01.2016 Сообщений: 663 Записей в блоге: 13
	24.03.2017, 20:56 [ТС]
	Я так понимаю, что во втором выражении в знаменателе нуль в нулевой степени при х=0? 0

@Cyborg Drone Модератор 10465 / 5761 / 3410 Регистрация: 17.08.2012 Сообщений: 17,515
	24.03.2017, 23:48
	Да. И, кроме того, при x=1 знаменатель обращается в 0, делить на 0 как-то нельзя, и функция при x=1 не определена (в данном случае ещё говорят, что значение функции стремится к бесконечности). Несмотря на всё это, определённый интеграл (правильнее сказать - определённый несобственный интеграл второго рода) сходится, иными словами, вопреки всему "очевидному", площадь между кривой и осью x конечна, и равна выше написано чему. 1

@SW Developer 97 / 93 / 81 Регистрация: 10.01.2016 Сообщений: 663 Записей в блоге: 13
	25.03.2017, 01:36 [ТС]
	И как быть, если в задании отрезок интегрирования [0, 7], а используя метод левых прямоугольников мы первым делом определяем Y при x=0? Добавлено через 3 минуты С параметром а проблем не будет, если в основном блоке программы прописать условие а<> 0. Добавлено через 31 секунду При вводе этого параметра с клавиатуры. Добавлено через 16 минут Получается, используя метод правых прямоугольников мы получаем численное значение интеграла, а используя метод левых прямоугольников - выходим за пределы интегрирования. 0

@bormant Модератор 7818 / 4637 / 2837 Регистрация: 22.11.2013 Сообщений: 13,159 Записей в блоге: 1
	25.03.2017, 13:29
	При делении 1 на бесконечно малое -- так и есть. 1

@SW Developer 97 / 93 / 81 Регистрация: 10.01.2016 Сообщений: 663 Записей в блоге: 13
	26.03.2017, 16:08 [ТС]
	А алгоритм по методу квадратуры Гаусса–Лежандра есть у кого-нибудь? 0

@SW Developer 97 / 93 / 81 Регистрация: 10.01.2016 Сообщений: 663 Записей в блоге: 13
	29.03.2017, 12:45 [ТС]
	Формулу использовал вот эту: $https://www.cyberforum.ru/cgi-bin/latex.cgi?<br /> I\approx \frac{b-a}{2}\left(f\left(\frac{a+b}{2}-\frac{b-a}{2\sqrt{3}} \right)+f\left(\frac{a+b}{2}+\frac{b-a}{2\sqrt{3}} \right)\right)<br />$ 0