Форум программистов, компьютерный форум, киберфорум
С++ для начинающих
Форум программистов и сисадминов Киберфорум > Форум программистов > C/C++ > C++ > С++ для начинающих
Войти
Регистрация
Восстановить пароль
Карта форума Темы раздела Блоги Сообщество Поиск Заказать работу  
 
Рейтинг 4.85/47: Рейтинг темы: голосов - 47, средняя оценка - 4.85
0 / 0 / 0
Регистрация: 27.02.2016
Сообщений: 39
1

Метод Монте-Карло для двойного интеграла с заданной точностью

04.12.2019, 20:14. Показов 8452. Ответов 4
Метки нет (Все метки)
Author24 — интернет-сервис помощи студентам
Здравствуйте, нужно вычислить интеграл https://www.cyberforum.ru/cgi-bin/latex.cgi?I = \int\limits_{ - 1}^1 {\int\limits_{ - 1 + \left| x \right|}^{1 - \left| x \right|} {\left| {x \cdot y \cdot (x + y)} \right|} } dydx методом Монте Карло с точностью https://www.cyberforum.ru/cgi-bin/latex.cgi?\varepsilon = {10^{ - 7}}. Есть значение интеграла, полученное в Mathcad'е https://www.cyberforum.ru/cgi-bin/latex.cgi?I = {\text{0}}{\text{.0916615}}. Есть рабочая программа, но она вычисляет максимум с точностью https://www.cyberforum.ru/cgi-bin/latex.cgi?\varepsilon = {10^{ - 5}}. Если точность задать больше, то условие в while не выполняется, следовательно происходит зацикливание. Что можете посоветовать в данном случае?
C++
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
33
34
35
36
37
38
39
40
41
42

#include <iostream>
#include <cmath>
using namespace std;
 
int main()
{
    double eps = 1e-5;
    double x, y,a=-1.0,b=1.0,c=-1.0,d=1.0;
    double Integr;
    double Integr_0 = 0.0;  // начальное значение для сравнения
    double summ = 0.0;      // для подсчета суммы
    int n = 0;              // для подсчета количества точек, попавших в заданную область
    int N = 500;            // первоначальное количество точек для вычисления
 
    
do {
        for (int i = 0; i < N; i++)
        {
            x = b + (a - b)*rand() / RAND_MAX;
            y = d + (c - d)*rand() / RAND_MAX;
            // если наша точка попадает в необходимую область
            // то вычисляем и накапливаем сумму, а также количество точек, попавших в данную область
            if (abs(x) + abs(y) <= 1)
            {
                summ += abs(x*y*(x + y));
                n++;
            }
        }
 
        // вычисляем наш интеграл. 2 - это площадь интегрирования. 
        Integr = 2 * summ / n;
        //меняем местами для сохранения в дальнейших вычислениях
        swap(Integr, Integr_0);
        // увеличиваем количество точек в 2 раза
        N *= 2;     
    } while (abs(Integr - Integr_0) > eps);
    // так как был обмен, то последний результат хранится в Integr_0
    cout << "n " << n << endl;
    cout << "N " << N << endl;
    cout <<"Integral "<< Integr_0 << endl;
    
}
0
Programming
Эксперт
94731 / 64177 / 26122
Регистрация: 12.04.2006
Сообщений: 116,782
 04.12.2019, 20:14
Ответы с готовыми решениями:

Программа для вычисления интеграла с заданной точностью методом Симпсона и методом Монте-Карло
Здравствуйте, подскажите как решить эту задачу, пожалуйста: Разработать программу для вычисления...

Метод Монте-Карло для вычисления площади фигуры
Суть такова, вычислить площадь фигуры, ограниченной кривыми x*y=a^2, x+y=5*a/2, где a вводится...

Вычисление интеграла методом Монте-Карло
Здравствуйте. Помогите в написании программы. Вычислить приближенно интеграл методом...

Вычисление интеграла методом Монте-Карло
Нужно вычислить интеграл, как показано на первой картинке внизу темы. На второй картинке сам...

4
случайный прохожий
3031 / 2062 / 626
Регистрация: 20.07.2013
Сообщений: 5,549
05.12.2019, 01:40 2
Попробуй использовать fabs вместо abs.
0
0 / 0 / 0
Регистрация: 27.02.2016
Сообщений: 39
05.12.2019, 07:16  [ТС] 3
fabs не помогло
0
3 / 2 / 1
Регистрация: 02.12.2018
Сообщений: 20
05.12.2019, 09:02 4
* Замените тип int на long long или unsigned long long. int очень быстро переполняется и уходит в минус. От этого и зацикливается.
* Вычисление с точностью 1e-6 измерсяеться минутами, а с точностью 1e-7 - часами...
* Точное решение интеграла быстрее всего 11/120 или 0.0916(6): вероятное решение

Вот код:
C++
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
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77

#include <iostream>
#include <cmath>
#include <iomanip>
#include <cstdint>
#include <chrono>
using namespace std;
using namespace chrono;
 
int main()
{
    constexpr double eps = 1e-5;
    constexpr double a = -1.0, b = 1.0;
    constexpr double c = -1.0, d = 1.0;
 
    double x, y;
 
    double Integr = 0.0;
    double Integr_0 = 0.0;  // начальное значение для сравнения
    double summ = 0.0;      // для подсчета суммы
    double distance;
 
    uint64_t n = 0;              // для подсчета количества точек, попавших в заданную область
    uint64_t N = 500;            // первоначальное количество точек для вычисления
 
 
    std::cout.precision(8);
    std::cout << std::fixed;
 
    cout << "              n                 N   |    Integr -   Integr_0| = distance" << endl;
 
    const time_point<system_clock> begin(chrono::system_clock::now());
 
    do {
        for (uint64_t i = 0; i < N; i++)
        {
            x = b + (a - b)*rand() / RAND_MAX;
            y = d + (c - d)*rand() / RAND_MAX;
 
            // если наша точка попадает в необходимую область
            // то вычисляем и накапливаем сумму, а также количество точек, попавших в данную область
            if (abs(x) + abs(y) <= 1)
            {
                summ += abs(x*y*(x + y));
                n++;
            }
        }
 
        // сохраняем предыдущее значение
        Integr_0 = Integr;
 
        // вычисляем наш интеграл. 2 - это площадь интегрирования.
        Integr = 2 * summ / n;
        distance = abs(Integr - Integr_0);
 
        std::cout << setw(15) << n << "   " << setw(15) << N << "   "
                  << '|' << setw(10) << Integr << " - " << setw(10) << Integr_0 << '|'
                  << " = " << setw(10) << distance << std::endl;
 
        // увеличиваем количество точек в 2 раза
        N *= 2;
 
    } while (distance > eps);
 
    const duration delta_time(system_clock::now() - begin);
 
    cout << "n " << n << endl;
    cout << "N " << N << endl;
    cout <<"Integral "<< Integr << endl;
 
 
   milliseconds milli = duration_cast< milliseconds>(delta_time) % seconds(1);
        seconds    ss = duration_cast<      seconds>(delta_time) % minutes(1);
        minutes    mm = duration_cast<      minutes>(delta_time) % hours(1);
          hours    hh = duration_cast<        hours>(delta_time);
    cout << "time spend " << hh.count() << ":" << mm.count() << ":"  << ss.count()
         << '.' << milli.count() << std::endl;
}
Результат при eps = 1e-5:
Кликните здесь для просмотра всего текста
Код
              n                 N   |    Integr -   Integr_0| = distance
            262               500   |0.09121678 - 0.00000000| = 0.09121678
            751              1000   |0.08893915 - 0.09121678| = 0.00227764
           1750              2000   |0.08756573 - 0.08893915| = 0.00137342
           3795              4000   |0.08712934 - 0.08756573| = 0.00043639
           7743              8000   |0.08979841 - 0.08712934| = 0.00266907
          15708             16000   |0.09135629 - 0.08979841| = 0.00155788
          31700             32000   |0.09119723 - 0.09135629| = 0.00015906
          63869             64000   |0.09117357 - 0.09119723| = 0.00002366
         128063            128000   |0.09153407 - 0.09117357| = 0.00036050
         255945            256000   |0.09146588 - 0.09153407| = 0.00006819
         512201            512000   |0.09154776 - 0.09146588| = 0.00008188
        1023550           1024000   |0.09151686 - 0.09154776| = 0.00003090
        2046898           2048000   |0.09159653 - 0.09151686| = 0.00007966
        4094697           4096000   |0.09160642 - 0.09159653| = 0.00000989
n 4094697
N 8192000
Integral 0.09160642
time spend 0:0:0.593


Результат при eps = 1e-6:
Кликните здесь для просмотра всего текста
Код
              n                 N   |    Integr -   Integr_0| = distance
            262               500   |0.09121678 - 0.00000000| = 0.09121678
            751              1000   |0.08893915 - 0.09121678| = 0.00227764
           1750              2000   |0.08756573 - 0.08893915| = 0.00137342
           3795              4000   |0.08712934 - 0.08756573| = 0.00043639
           7743              8000   |0.08979841 - 0.08712934| = 0.00266907
          15708             16000   |0.09135629 - 0.08979841| = 0.00155788
          31700             32000   |0.09119723 - 0.09135629| = 0.00015906
          63869             64000   |0.09117357 - 0.09119723| = 0.00002366
         128063            128000   |0.09153407 - 0.09117357| = 0.00036050
         255945            256000   |0.09146588 - 0.09153407| = 0.00006819
         512201            512000   |0.09154776 - 0.09146588| = 0.00008188
        1023550           1024000   |0.09151686 - 0.09154776| = 0.00003090
        2046898           2048000   |0.09159653 - 0.09151686| = 0.00007966
        4094697           4096000   |0.09160642 - 0.09159653| = 0.00000989
        8192265           8192000   |0.09159687 - 0.09160642| = 0.00000955
       16389017          16384000   |0.09163554 - 0.09159687| = 0.00003867
       32769576          32768000   |0.09165066 - 0.09163554| = 0.00001512
       65535286          65536000   |0.09164284 - 0.09165066| = 0.00000783
      131078815         131072000   |0.09165979 - 0.09164284| = 0.00001695
      262140246         262144000   |0.09166859 - 0.09165979| = 0.00000880
      524281599         524288000   |0.09166323 - 0.09166859| = 0.00000536
     1048545399        1048576000   |0.09166440 - 0.09166323| = 0.00000117
     2097077721        2097152000   |0.09166610 - 0.09166440| = 0.00000170
     4194273976        4194304000   |0.09166681 - 0.09166610| = 0.00000071
n 4194273976
N 8388608000
Integral 0.09166681
time spend 0:8:54.365


eps = 1e-7 (всё, что есть):
Кликните здесь для просмотра всего текста
Код
              n                 N   |    Integr -   Integr_0| = distance
            262               500   |0.09121678 - 0.00000000| = 0.09121678
            751              1000   |0.08893915 - 0.09121678| = 0.00227764
           1750              2000   |0.08756573 - 0.08893915| = 0.00137342
           3795              4000   |0.08712934 - 0.08756573| = 0.00043639
           7743              8000   |0.08979841 - 0.08712934| = 0.00266907
          15708             16000   |0.09135629 - 0.08979841| = 0.00155788
          31700             32000   |0.09119723 - 0.09135629| = 0.00015906
          63869             64000   |0.09117357 - 0.09119723| = 0.00002366
         128063            128000   |0.09153407 - 0.09117357| = 0.00036050
         255945            256000   |0.09146588 - 0.09153407| = 0.00006819
         512201            512000   |0.09154776 - 0.09146588| = 0.00008188
        1023550           1024000   |0.09151686 - 0.09154776| = 0.00003090
        2046898           2048000   |0.09159653 - 0.09151686| = 0.00007966
        4094697           4096000   |0.09160642 - 0.09159653| = 0.00000989
        8192265           8192000   |0.09159687 - 0.09160642| = 0.00000955
       16389017          16384000   |0.09163554 - 0.09159687| = 0.00003867
       32769576          32768000   |0.09165066 - 0.09163554| = 0.00001512
       65535286          65536000   |0.09164284 - 0.09165066| = 0.00000783
      131078815         131072000   |0.09165979 - 0.09164284| = 0.00001695
      262140246         262144000   |0.09166859 - 0.09165979| = 0.00000880
      524281599         524288000   |0.09166323 - 0.09166859| = 0.00000536
     1048545399        1048576000   |0.09166440 - 0.09166323| = 0.00000117
     2097077721        2097152000   |0.09166610 - 0.09166440| = 0.00000170
     4194273976        4194304000   |0.09166681 - 0.09166610| = 0.00000071
     8388573769        8388608000   |0.09166597 - 0.09166681| = 0.00000084
    16777343125       16777216000   |0.09166656 - 0.09166597| = 0.00000059
    33554619698       33554432000   |0.09166674 - 0.09166656| = 0.00000018
    67108922918       67108864000   |0.09166640 - 0.09166674| = 0.00000034
...
0
693 / 303 / 99
Регистрация: 04.07.2014
Сообщений: 846
05.12.2019, 14:16 5
Несколько важных замечаний:

1. Генератор случайных чисел

rand обладает очень маленьким периодом, и для генерации миллиардов случайных чисел не подходит. Для Метода Монте-Карло очень важен правильный выбор генератора псевдослучайных чисел.

Так что смотрим в сторону https://ru.cppreference.com/w/cpp/numeric/random

2. Суммирование большого числа чисел типа double

Складывать несколько миллионов чисел типа double просто так нельзя.

Смотрим:

3. Точность метода и распределение случайной величины

https://habr.com/ru/post/274975/ - крайне полезная статья.

Итог

Для достижения требуемой точности надо объединить все три решения.

З.Ы.: Вечером попробую что-нибудь сделать
0
05.12.2019, 14:16
IT_Exp
Эксперт
87844 / 49110 / 22898
Регистрация: 17.06.2006
Сообщений: 92,604
 05.12.2019, 14:16
Помогаю со студенческими работами здесь

Вычисление интеграла методом Монте-Карло
Здравствуйте! Задача такая: пользователь в текстовом файле задает координаты точек (любое...

Метод Монте-Карло для решение систем линейных алгебраических уравнений
Здравствуйте форумчане, у меня, наверное, проблема больше с математикой. Попытался написать...

Решения кратного интеграла методом Монте Карло на С++
Помогите пожалуста решить тройной интеграл методом Монте Карло.... Нужно написать програму на С ...

Вычисление интеграла геометрическим методом Монте-Карло
Всем доброго времени суток. В универе дали задание: вычислить интеграл...


Искать еще темы с ответами

Или воспользуйтесь поиском по форуму:
5
Ответ Создать тему
КиберФорум - форум программистов, компьютерный форум, программирование
Поддержать - Реклама - Обратная связь
Powered by vBulletin
Copyright ©2000 - 2024, CyberForum.ru